simple geometric model

1. description of all variables

labeled_cartoon.png
\begin{equation} F=\sigma\left(S_{0}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}
description value units reference type python variable associated equation LaTeX variable
invagination height   \(nm\)   input (variable) Hi   \(H_{i}\)
membrane tension 0.001 \(pN/nm\) hassinger_design_2017 fixed sigma   \(\sigma\)
total plane radius 2000 \(nm\)   fixed At   \(A_{t}\)
final vesicle radius 50 \(nm\) kaksonen_mechanisms_2018 fixed vesicleradius   \(R_{v}\)
coat bending rigidity 160 \(pN nm\) dmitrieff_membrane_2015 fixed K   \(K\)
spontaneous curvature 1/50 \(nm^{-1}\) hassinger_design_2017 fixed Ci   \(C_{i}\)
droplet-cytosol surface tension 0.005 \(pN/nm\) jawerth_salt-dependent_2018 fixed gammad   \(\gamma_{d}\)
droplet-membrane wetting tension 0.002 \(pN/nm\)   fixed gammam   \(\gamma_{m}\)
initial droplet plane radius       fixed      
droplet-invagination surface tension 0.002 \(pN/nm\)   fixed gammai   \(\gamma_{i}\)
turgor pressure 0.1 \(pN/nm^2\) ma_mechanics_2019 fixed P   \(P\)
relative invagination height   unitless   input (variable) h rel-h \(h\)
droplet volume   \(nm^3\)   intermediate output Vd   \(V_{d}\)
invagination surface area   \(nm^2\)   intermediate output Si invagination-area \(S_{i}\)
initial droplet spherical radius 600 \(nm\)   intermediate output dropradius    
droplet-flat membrane contact area   \(nm^2\)   intermediate output Sm md-area \(S_{m}\)
invagination spherical radius   \(nm\)   intermediate output Ri invagination-radius \(R_{i}\)
droplet-cytosol contact area   \(nm^2\)   intermediate output Sd dc-contact-area \(S_{d}\)
droplet contact angle   rad   intermediate output thetad droplet-contact-angle \(\theta_{d}\)
invagination contact angle   rad   intermediate output thetai invagination-angle \(\theta_{i}\)
invagination plane radius   \(nm\)   intermediate output Ai invagination-plane-radius \(A_{i}\)
droplet + invagination plane radius   \(nm\)   intermediate output Ad droplet-plane-radius \(A_{d}\)
outside membrane surface area   \(nm^2\)   intermediate output So outside-membrane-area \(S_{o}\)
invagination volume   \(nm^3\)   intermediate output Vi invagination-volume \(V_{i}\)
droplet + invagination volume   \(nm^3\)   intermediate output Vb di-volume \(V_{b}\)
droplet + invagination height   \(nm\)   intermediate output Hd di-height \(H_{d}\)
droplet + invagination spherical radius   \(nm\)   intermediate output Rd di-cap-radius \(R_{d}\)
free energy   \(pN \cdot nm\)   output F f-forces f-shapes \(F\)
membrane tension energy   \(pN \cdot nm\)   output membranetension    
bending energy   \(pN \cdot nm\)   output bendingenergy    
surface tension energy   \(pN \cdot nm\)   output surfacetension    
wetting energy   \(pN \cdot nm\)   output wettingenergy    
turgor pressure energy   \(pN \cdot nm\)   output turgorpressure    
membrane curvature   \(nm^{-1}\)   output 2÷Ri    

4.114 \(pN \cdot nm\) = 1 \(k_{B}T\) noauthor_kt_2021

2. equations

2.1. summary

labeled_cartoon.png

arranged by force contributions:

F = membrane tension + bending energy + droplet-cytosol surface tension -
droplet-membrane surface tension - droplet-invagination surface tension + turgor
pressure 
\begin{equation} F=\sigma\left(S_{o}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}

rearranged by shapes:

\begin{equation} F= \sigma S_{o} + \left(\sigma - \gamma_{m} \right) S_{m} + \gamma_{d} S_{d} + \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} + P V_{i} + \left(\sigma - \gamma_{i} \right) S_{i} \end{equation}

broken down to shape functions and constant \(C\):

\begin{equation} F(\Delta S_{o}, \Delta S_{m}, S_{d}, R_{i}, V_{i})= F_{1}(\Delta S_{o}) + F_{2}(\Delta S_{m}) + F_{3}(S_{d}) + F_{4}(R_{i}) + F_{5}(V_{i}) + C \end{equation} \begin{equation} F_{1}(\Delta S_{o}) = \sigma \Delta S_{o} \end{equation} \begin{equation} \Delta S_{o}=\alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \label{alpha} \alpha = - \left(\frac{\pi\left(\frac{\gamma_{m}}{\gamma_{d}}+1\right)^{3}} {\left(-\frac{\gamma_{m}}{\gamma_{d}}+1\right)\left(\frac{\gamma_{m}}{\gamma_{d}}+2\right)^{2}}\right)^{1 / 3} \end{equation} \begin{equation} \label{} F_{2}(\Delta S_{m}) = \left(\sigma - \gamma_{m} \right) \Delta S_{m} \end{equation} \begin{equation} \Delta S_{m}= - \alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \label{} F_{3}(S_{d}) = \gamma_{d} S_{d} \end{equation} \begin{equation} S_{d}=\epsilon \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \epsilon = \left(\frac{8 \pi}{\left(-\frac{\gamma_{m}}{\gamma_{d}}+1\right) \left(\frac{\gamma_{m}}{\gamma_{d}}+2\right)^{2}}\right)^{1 / 3} \end{equation} \begin{equation} \label{} \begin{array}{l} F_{4}(R_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} \end{array} \end{equation} \begin{equation} R_{i} = \frac{\sqrt{S_i}}{2 \sqrt{\pi } h} \end{equation} \begin{equation} \label{} F_{4}(h) = 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^{2} \end{equation} \begin{equation} \label{} F_{5}(V_{i}) = P V_{i} \end{equation} \begin{equation} \label{} V_{i}=\eta \left(3 h-2 h^{3}\right) \end{equation} \begin{equation} \label{} \eta=\frac{S_i^{3/2}}{6 \sqrt{\pi }} \end{equation} \begin{equation} \label{} C = \sigma \pi A_t^2+S_i \left(\gamma _m-\gamma _i\right) \end{equation}

all combined into full free energy equation:

\begin{equation} F(h)=\left(\gamma_{m} \alpha+\gamma_{d} \epsilon\right)\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3}+\left(\sigma-\gamma_{m}\right) S_{i} h^{2} +2 K\left(C_{i} \sqrt{S}_{i}-2 \sqrt{\pi} h\right)^{2}+P \eta \left(3 h-2 h^{3}\right)+C \end{equation}

expanded polynomial:

\begin{equation} F(h)= \left(\gamma_{m} \alpha+\gamma_{d} \epsilon\right) \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} - 2 P \eta h^{3}+\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h^{2} +\left(3 P \eta-8 K C_{i} \sqrt{\pi S_{i}}\right) h+B \end{equation} \begin{equation} B=C+2 K C_{i}^{2} S_{i} \end{equation}

2.2. individual shape solutions

2.2.1. relative invagination height

h.jpg
\begin{equation} \label{rel-h} h = \frac{H_{i}}{\left(\frac{S_{i}}{\pi}\right)^{1/2}} \end{equation} \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation}

2.2.2. droplet-membrane contact angle

2020-09-21_19-17-56_screenshot.png
Figure 1: contact angle cartoon contact-angle
\begin{equation} \begin{array}{l} \gamma_{S G}=\gamma_{S L}+\gamma_{L G} \cos \theta \\ \sigma+\gamma_{m}=\sigma+\gamma_{d} \cos \theta_{d} \\ \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \\ \sin \theta_{d}=\sqrt{1 - \cos \theta_{d}^{2}} \\ \sin \theta_{d}=\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}} \end{array} \end{equation} \begin{equation} \label{droplet-contact-angle} \theta_{d}=\arccos \left(\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

2.2.3. outside membrane area

\begin{equation} \label{outside-membrane-area} S_{0}=\pi\left(A_{t}^{2}-A_{d}^{2}\right) \end{equation}

2.2.4. circular plane radii

\begin{equation} \begin{array}{l} \sin \theta=\frac{a}{r} \\ a=r \sin \theta} \\ A_{i}=R_{i} \sin \theta_{i} \end{array} \end{equation} \begin{equation} \label{droplet-plane-radius} A_{d}=R_{d} \sin \theta_{d} \end{equation}

2.2.5. droplet-cytosol contact area

\begin{equation} \label{dc-contact-area} \begin{array}{l} A=2 \pi r^{2}(1-\cos \theta) \\ S_{d}=2 \pi R_{d}^{2}\left(1-\cos \theta_{d}\right) \end{array} \end{equation}

2.2.6. invagination surface area

\begin{equation} \label{invagination-area} S_{i}=4 \pi R_{v}^{2} \end{equation}

2.2.7. invagination volume

\begin{equation} \label{invagination-volume-a-h} \begin{array}{l} V=\frac{1}{6} \pi h\left(3 a^{2}+h^{2}\right) \\ V_{i}=\frac{1}{6} \pi H_{i}\left(3 A_{i}^{2}+H_{i}^{2}\right) \end{array} \end{equation}

alternately:

\begin{equation} \label{invagination-volume} \begin{array}{l} V=\frac{\pi}{3} h^{2} \left(3 r - h \right) \\ V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{array} \end{equation}

2.2.8. invagination plane radius

\begin{equation} \begin{array}{l} A=\pi\left(a^{2}+h^{2}\right) \\ S_{i}=\pi\left(A_{i}^{2}+H_{i}^{2}\right) \\ A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{array} \end{equation} \begin{equation} \label{invagination-plane-radius} A_{i}=\sqrt{\frac{S_{i}}{\pi} \left(1 - h^{2}\right)} \end{equation}

2.2.9. invagination radius

\begin{equation} \label{invagination-radius} \begin{aligned} A &=2 \pi r h \\ S_{i} &=2 \pi R_{i} H_{i} \\ R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{aligned} \end{equation}

2.2.10. invagination contact angle

\begin{equation} \label{invagination-angle} \begin{array}{l} S_{i}=2 \pi R_{i}^{2}\left(1-\cos \theta_{i}\right) \\ \theta_{i}=\arccos \left(1-\frac{S_{i}}{2 \pi R_{i}^{2}}\right) \end{array} \end{equation}

2.2.11. TODO droplet volume

  • should depend on initial plane radius

    \begin{equation} \label{d-volume} V_{b}=V_{d}+V_{i} \end{equation}

2.2.12. droplet + invagination volume

\begin{equation} \label{di-volume} V_{b}=V_{d}+V_{i} \end{equation}

2.2.13. droplet + invagination spherical cap radius

\begin{equation} \label{di-cap-radius} \begin{array}{l} V=\frac{\pi}{3} r^{3}(2+\cos \theta)(1-\cos \theta)^{2} \\ V_{b}=\frac{\pi}{3} R_{d}^{3}\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2} \\ R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{array} \end{equation}

2.2.14. droplet + invagination height

\begin{equation} \label{di-height} \begin{array}{l} A=2\pi r h \\ S_{d}=2\pi R_{d} H_{d} \\ H_{d}=\frac{S_{d}}{2 \pi R_{d}} \end{array} \end{equation}

2.2.15. flat membrane-droplet contact area

\begin{equation} \label{md-area} S_{m}=\pi\left(A_{d}^{2}-A_{i}^{2}\right) \end{equation}

2.3. separate force components (original version)

\begin{equation} \label{f-forces} F=\sigma\left(S_{o}+S_{m}+S_{i}\right)+\frac{K}{2} S_{i}\left(\frac{2}{R_{i}}-2C_{i}\right)^{2}+\gamma_{d} S_{d}-\gamma_{m} S_{m}-\gamma_{i} S_{i}+P V_{i} \end{equation}

2.4. fully substituted

\begin{equation} \label{free-energy} F=\sigma(\pi(A_{t}^{2}-(((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2})+(\pi((((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2}-(\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}))+(4 \pi R_{v}^{2}))+ \frac{K}{2} (4 \pi R_{v}^{2})(\frac{2}{(\frac{(4 \pi R_{v}^{2})}{2 \pi \left(\frac{S_{i} h}{\pi}\right)})}-2C_{i})^{2}+ \gamma_{d} (2 \pi R_{d}^{2}(1-(\frac{\gamma_{m}}{\gamma_{d}})))- \gamma_{m} (\pi((((\frac{3 (V_{d}+(\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})))}{\pi(2+(\frac{\gamma_{m}}{\gamma_{d}}))(1-(\frac{\gamma_{m}}{\gamma_{d}}))^{2}})^{1 / 3}) (\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}))^{2}-(\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}))- \gamma_{i} (4 \pi R_{v}^{2})+ P (\frac{1}{6} \pi \left(\frac{S_{i} h}{\pi}\right)(3 (\sqrt{\frac{(4 \pi R_{v}^{2})}{\pi}-\left(\frac{S_{i} h}{\pi}\right)^{2}})^{2}+\left(\frac{S_{i} h}{\pi}\right)^{2})) \end{equation}

2.5. separate shape force components

\begin{equation} F(\Delta S_{o}, \Delta S_{m}, S_{d}, R_{i}, V_{i})= F_{1}(\Delta S_{o}) + F_{2}(\Delta S_{m}) + F_{3}(S_{d}) + F_{4}(R_{i}) + F_{5}(V_{i}) + C \end{equation} \begin{equation} \label{f-shapes} F= \sigma S_{o} + \left(\sigma - \gamma_{m} \right) S_{m} + \gamma_{d} S_{d} + \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} + P V_{i} + D \end{equation}

2.5.1. full derivations

note: there may be some residual messyness since I caught a mistake in the sign of \(\gamma_{d}\) in droplet-contact-angle and corrected it after writing all these equations

  1. \(F_{1}(\Delta S_{o})\)
    \begin{equation} F_{6}(S_{o}) = \sigma S_{o} \end{equation}

    substituting outside-membrane-area:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-A_{d}^{2}\right) \end{equation}

    substituting droplet-plane-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-(R_{d} \sin \theta_{d})^{2}\right) \end{equation}

    substituting di-cap-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}\right) \end{equation}

    substituting droplet-contact-angle:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 V_{b}}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting di-volume:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 (V_{d}+V_{i})}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting invagination-volume:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    substituting invagination-radius:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi\left(A_{t}^{2}-\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}\right) \end{equation}

    "simplified" in mathematica:

    \begin{equation} F_{6}(S_{o}) = \sigma \pi \left(A_t^2+\frac{\left(\gamma _m^2-\gamma _d^2\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(\frac{\gamma _m}{\gamma _d}+2\right) \left(\gamma _d-\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d^2}\right) \end{equation}

    rearranged to isolate constants vs. functions of \(H_{i}\):

    endocytosis phase separation V2 - Page 5.png
    \begin{equation} F_{6}(S_{o}) = \sigma S_{o} \end{equation} \begin{equation} S_{o} = \Delta S_{o} + \pi A_{t}^{2} \end{equation} \begin{equation} F_{6}(S_{o}) = F_{1}(\Delta S_{o}) + \sigma \pi A_{t}^{2} \end{equation} \begin{equation} F_{1}(\Delta S_{o}) = \sigma \Delta S_{o} \end{equation} \begin{equation} \Delta S_{o} = S_{o} - \pi A_{t}^{2} \end{equation} \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi H_{i}^{3}+\frac{3}{2} S_{i} H_{i}+3 V_{d}\right)^{2 / 3} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{3}+\frac{3}{2} S_{i} h \left(\frac{S_{i}}{\pi}\right)^{1/2}+3 V_{d}\right)^{2 / 3} \end{equation}

    simplified:

    \begin{equation} \Delta S_{o}=\left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \Delta S_{o}=\alpha \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \alpha = \left(\frac{\pi\left(\gamma_{m}^{2}-\gamma_{d}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \end{equation}
  2. \(F_{2}(\Delta S_{m})\)

    gathering from free-energy:

    \begin{equation} \label{} F_{7}(S_{m}) = \left(\sigma - \gamma_{m} \right) S_{m} \end{equation}

    from md-area:

    \begin{equation} \label{md-area} S_{m}=\pi\left(A_{d}^{2}-A_{i}^{2}\right) \end{equation}

    substitute droplet-plane-radius:

    \begin{equation} A_{d}=R_{d} \sin \theta_{d} \end{equation}

    and invagination-plane-radius:

    \begin{equation} A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(R_{d} \sin \theta_{d}\right)^{2}-\frac{S_{i}}{\pi}+H_{i}^{2}\right) \end{equation}

    rearranged:

    \begin{equation} S_{m}=\pi\left(\left(R_{d} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute di-cap-radius:

    \begin{equation} R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute di-volume:

    \begin{equation} V_{b}=V_{d}+V_{i} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+V_{i} \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute invagination-volume:

    \begin{equation} V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute invagination-radius:

    \begin{equation} R_{i} =\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right) \right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \sin \theta_{d}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    substitute droplet-contact-angle:

    \begin{equation} \begin{array}{l} \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \\ \sin \theta_{d}=\sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}} \end{array} \end{equation}

    to get:

    \begin{equation} S_{m}=\pi\left(\left(\left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right) \right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{1 / 3} \sqrt{1 - (\frac{\gamma_{m}}{\gamma_{d}})^{2}}\right)^{2}+H_{i}^{2}\right) - S_{i} \end{equation}

    simplified in mathematica:

    \begin{equation} S_{m}=\pi \left(\frac{\left(\gamma _d^2-\gamma _m^2\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(2-\frac{\gamma _m}{\gamma _d}\right) \left(\gamma _d+\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d^2}+H_i^2\right)-S_i \end{equation}

    manually rearranged:

    endocytosis phase separation V2 - Page 6.png
    \begin{equation} \label{} F_{7}(S_{m}) = \left(\sigma - \gamma_{m} \right) S_{m} \end{equation} \begin{equation} S_{m}= \Delta S_{m} - S_{i} \end{equation} \begin{equation} \label{} F_{7}(S_{m}) = F_{2}(\Delta S_{m}) - \left(\sigma - \gamma_{m} \right) S_{i} \end{equation} \begin{equation} \label{} F_{2}(\Delta S_{m}) = \left(\sigma - \gamma_{m} \right) \Delta S_{m} \end{equation} \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi H_{i}^{3}+\frac{3}{2} S_{i} H_{i}+3 V_{d}\right)^{2 / 3}+\pi H_{i}^{2} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3}\left(-\pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{3}+\frac{3}{2} S_{i} \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)+3 V_{d}\right)^{2 / 3}+ \pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^{2} \end{equation}

    simplified:

    \begin{equation} \Delta S_{m}=\left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \Delta S_{m}= \beta \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} + S_{i} h^{2} \end{equation} \begin{equation} \beta = \left(\frac{\pi\left(\gamma_{d}^{2}-\gamma_{m}^{2}\right)^{3}}{\left(2 \gamma_{d}+\gamma_{m}\right)^{2}\left(\gamma_{d}-\gamma_{m}\right)^{4}}\right)^{1 / 3} \end{equation}
  3. \(F_{3}(S_{d})\)

    gathering from free-energy:

    \begin{equation} \label{} F_{3}(S_{d}) = \gamma_{d} S_{d} \end{equation}

    from dc-contact-area:

    \begin{equation} S_{d}=2 \pi R_{d}^{2}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute di-cap-radius

    \begin{equation} R_{d}=\left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{1 / 3} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 V_{b}}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute di-volume

    \begin{equation} V_{b}=V_{d}+V_{i} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+V_{i}\right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute invagination-volume

    \begin{equation} V_{i}=\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right) \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\cos \theta_{d}\right)\left(1-\cos \theta_{d}\right)^{2}}\right)^{2 / 3}\left(1-\cos \theta_{d}\right) \end{equation}

    substitute droplet-contact-angle

    \begin{equation} \cos \theta_{d}=\frac{\gamma_{m}}{\gamma_{d}} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 R_{i} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{2 / 3}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

    substitute invagination-radius

    \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    to get:

    \begin{equation} S_{d}=2 \pi \left(\frac{3 \left(V_{d}+\frac{\pi}{3} H_{i}^{2} \left(3 \frac{S_{i}}{2 \pi H_{i}} - H_{i} \right)\right)}{\pi\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}}\right)^{2 / 3}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right) \end{equation}

    simplified in mathematica:

    \begin{equation} S_{d}= \frac{2 \pi \left(\gamma _d-\gamma _m\right) \left(\frac{\gamma _d^2 \left(3 V_d+\pi H_i^2\left(\frac{3 S_i}{2 \pi H_i}-H_i\right)\right)}{\pi \left(\frac{\gamma _m}{\gamma _d}+2\right) \left(\gamma _d-\gamma _m\right)^{2}}\right)^{2/3}}{\gamma _d} \end{equation}

    manually rearranged:

    20210104_sd.png
    \begin{equation} S_{d}=\left(\frac{8 \pi}{\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}\left(1-\frac{\gamma_{m}}{\gamma_{d}}\right)}\right)^{1 / 3}\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} S_{d}=\epsilon \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation} \begin{equation} \epsilon = \left(\frac{8 \pi}{\left(2+\frac{\gamma_{m}}{\gamma_{d}}\right)^{2}\left(1-\frac{\gamma_{d}}{\gamma_{m}}\right)}\right)^{1 / 3} \end{equation}
  4. \(F_{4}(R_{i})\)

    gathering from f-forces:

    \begin{equation} \label{} \begin{array}{l} F_{4}(R_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{R_{i}} - 2 C_{i}\right)^{2} \end{array} \end{equation}

    from invagination-radius

    \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} R_{i} = \frac{S_{i}}{2 \pi h \left(\frac{S_{i}}{\pi}\right)^{1/2}} \end{equation}

    simplified:

    \begin{equation} R_{i} = \frac{\sqrt{S_i}}{2 \sqrt{\pi } h} \end{equation}
    1. \(F_{4}(h)\)

      substitute invagination-radius

      \begin{equation} R_{i} &=\frac{S_{i}}{2 \pi H_{i}} \end{equation}

      to get:

      \begin{equation} \label{} \begin{array}{l} F_{4}(H_{i}) = \frac{K}{2} S_{i} \left(\frac{2}{\frac{S_{i}}{2 \pi H_{i}}} - 2 C_{i}\right)^{2} \end{array} \end{equation}

      simplified:

      \begin{equation} \label{} F_{4}(H_{i}) = \frac{2 K}{S_i} \left(C_i S_i-2 \pi H_i\right)^{2} \end{equation}

      replacing invagination height from rel-h:

      \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \label{} F_{4}(h) = \frac{2 K}{S_i} \left(C_i S_i-2 \pi \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)\right)^{2} \end{equation}

      simplified:

      \begin{equation} \label{} F_{4}(h) = 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^2 \end{equation}
  5. \(F_{5}(V_{i})\)

    gathering from f-forces:

    \begin{equation} \label{} F_{5}(V_{i}) = P V_{i} \end{equation}

    from invagination-volume-a-h:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(3 A_{i}^{2}+H_{i}^{2}\right) \end{array} \end{equation}

    substitute invagination-plane-radius

    \begin{equation} \label{} \begin{array}{l} A_{i}=\sqrt{\frac{S_{i}}{\pi}-H_{i}^{2}} \end{array} \end{equation}

    to get:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(3 \left(\frac{S_{i}}{\pi}-H_{i}^{2}\right)+H_{i}^{2}\right) \end{array} \end{equation}

    simplified:

    \begin{equation} \label{} \begin{array}{l} V_{i}=\frac{1}{6} \pi H_{i}\left(\frac{3 S_i}{\pi }-2 H_{i}^2\right) \end{array} \end{equation}

    replacing invagination height from rel-h:

    \begin{equation} H_{i} = h \left(\frac{S_{i}}{\pi}\right)^{1/2} \end{equation} \begin{equation} \label{} V_{i}=\frac{1}{6} \pi h \left(\frac{S_{i}}{\pi}\right)^{1/2}\left(\frac{3 S_i}{\pi }-2 \left(h \left(\frac{S_{i}}{\pi}\right)^{1/2}\right)^2\right) \end{equation}

    simplified:

    \begin{equation} \label{} V_{i}= \frac{h \left(3-2 h^2\right) S_i^{3/2}}{6 \sqrt{\pi }} \end{equation} \begin{equation} \label{} V_{i}=\eta \left(3 h-2 h^{3}\right) \end{equation} \begin{equation} \label{} \eta=\frac{S_i^{3/2}}{6 \sqrt{\pi }} \end{equation}
  6. \(D\): everything else that doesn't change with \(H_{i}\)
    \begin{equation} \label{D} D = \sigma S_{i} - \gamma_{i} S_{i} \end{equation}
  7. \(C\): additional constants grouped out from \(\Delta\) terms
    \begin{equation} \label{} C = D + \sigma \pi A_{t}^{2} - \left(\sigma - \gamma_{m} \right) S_{i} \end{equation} \begin{equation} \label{} D = \sigma S_{i} - \gamma_{i} S_{i} \end{equation}

    simplified:

    \begin{equation} \label{} C = \sigma \pi A_t^2+S_i \left(\gamma _m-\gamma _i\right) \end{equation}

2.5.2. DONE substituted after isolating \(h\)

  • could anything be factored out more first?
\begin{equation} \label{} F(h) = \sigma \alpha \xi + \left(\sigma - \gamma_{m} \right) \left( \beta \xi + S_{i} h^{2}\right) + \gamma_{d} \epsilon \xi + 2 K \left(C_i \sqrt{S_i}-2 \sqrt{\pi } h\right)^{2} + P \eta \left(3 h-2 h^{3}\right) + C \end{equation} \begin{equation} \label{} \xi = \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} \end{equation}
20201203_fh (1).png
\begin{equation} F(h)=\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right)\left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3}+\left(\sigma-\gamma_{m}\right) S_{i} h^{2} +2 K\left(C_{i} \sqrt{S}_{i}-2 \sqrt{\pi} h\right)^{2}+P \eta \left(3 h-2 h^{3}\right)+C \end{equation} \begin{equation} F(h)= \left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) \left(\frac{S_{i}^{3 / 2}}{\sqrt{\pi}}\left(\frac{3}{2} h-h^{3}\right)+3 V_{d}\right)^{2/3} - 2 P \eta h^{3}+\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h^{2} +\left(3 P \eta-8 C_{i} \sqrt{\pi S_{i}}\right) h+B \end{equation} \begin{equation} B=C+2 K C_{i}^{2} S_{i} \end{equation}

2.5.3. testing sign of alpha terms

import numpy as np
import math
def cube_root(i):
    return math.pow(abs(i),float(1)/3) * (1,-1)[i<0]

def alpha_so(gamma_m, gamma_d):
    alpha = cube_root((np.pi*(gamma_m**2 - gamma_d**2)**3)/((2*gamma_d+gamma_m)**2 * (gamma_d-gamma_m)**4))
    return alpha

for gamma_m in np.linspace(0, 200, 9):
    for gamma_d in np.linspace(0, 200, 9):
        print('gamma_m',gamma_m,'gamma_d',gamma_d,alpha_so(gamma_m, gamma_d))
gamma_m 0.0 gamma_d 0.0 nan
gamma_m 0.0 gamma_d 25.0 -0.9226350743220142
gamma_m 0.0 gamma_d 50.0 -0.9226350743220142
gamma_m 0.0 gamma_d 75.0 -0.9226350743220142
gamma_m 0.0 gamma_d 100.0 -0.9226350743220142
gamma_m 0.0 gamma_d 125.0 -0.9226350743220142
gamma_m 0.0 gamma_d 150.0 -0.9226350743220142
gamma_m 0.0 gamma_d 175.0 -0.9226350743220142
gamma_m 0.0 gamma_d 200.0 -0.9226350743220142
gamma_m 25.0 gamma_d 0.0 1.4645918875615231
gamma_m 25.0 gamma_d 25.0 nan
gamma_m 25.0 gamma_d 50.0 -1.5026501396568157
gamma_m 25.0 gamma_d 75.0 -1.2706753068765375
gamma_m 25.0 gamma_d 100.0 -1.1735039002840684
gamma_m 25.0 gamma_d 125.0 -1.1192302015380546
gamma_m 25.0 gamma_d 150.0 -1.084415611760493
gamma_m 25.0 gamma_d 175.0 -1.0601370721917294
gamma_m 25.0 gamma_d 200.0 -1.042222647325497
gamma_m 50.0 gamma_d 0.0 1.4645918875615231
gamma_m 50.0 gamma_d 25.0 1.7436710272644398
gamma_m 50.0 gamma_d 50.0 nan
gamma_m 50.0 gamma_d 75.0 -1.830739859451904
gamma_m 50.0 gamma_d 100.0 -1.5026501396568157
gamma_m 50.0 gamma_d 125.0 -1.356188576761231
gamma_m 50.0 gamma_d 150.0 -1.2706753068765375
gamma_m 50.0 gamma_d 175.0 -1.214010595446712
gamma_m 50.0 gamma_d 200.0 -1.1735039002840684
gamma_m 75.0 gamma_d 0.0 1.4645918875615231
gamma_m 75.0 gamma_d 25.0 1.590205608287594
gamma_m 75.0 gamma_d 50.0 2.001188208394222
gamma_m 75.0 gamma_d 75.0 nan
gamma_m 75.0 gamma_d 100.0 -2.0727783992021025
gamma_m 75.0 gamma_d 125.0 -1.6820324801561286
gamma_m 75.0 gamma_d 150.0 -1.5026501396568157
gamma_m 75.0 gamma_d 175.0 -1.3955026949998752
gamma_m 75.0 gamma_d 200.0 -1.3231738439535343
gamma_m 100.0 gamma_d 0.0 1.4645918875615231
gamma_m 100.0 gamma_d 25.0 1.5377251238700236
gamma_m 100.0 gamma_d 50.0 1.7436710272644398
gamma_m 100.0 gamma_d 75.0 2.208757298511275
gamma_m 100.0 gamma_d 100.0 nan
gamma_m 100.0 gamma_d 125.0 -2.2692052337015594
gamma_m 100.0 gamma_d 150.0 -1.830739859451904
gamma_m 100.0 gamma_d 175.0 -1.6263744927117951
gamma_m 100.0 gamma_d 200.0 -1.5026501396568157
gamma_m 125.0 gamma_d 0.0 1.4645918875615231
gamma_m 125.0 gamma_d 25.0 1.512803489134373
gamma_m 125.0 gamma_d 50.0 1.6429054603976958
gamma_m 125.0 gamma_d 75.0 1.8801889207945017
gamma_m 125.0 gamma_d 100.0 2.384131644400035
gamma_m 125.0 gamma_d 125.0 nan
gamma_m 125.0 gamma_d 150.0 -2.436744690673985
gamma_m 125.0 gamma_d 175.0 -1.959079850097313
gamma_m 125.0 gamma_d 200.0 -1.7343631139416589
gamma_m 150.0 gamma_d 0.0 1.4645918875615231
gamma_m 150.0 gamma_d 25.0 1.498872430465395
gamma_m 150.0 gamma_d 50.0 1.590205608287594
gamma_m 150.0 gamma_d 75.0 1.7436710272644396
gamma_m 150.0 gamma_d 100.0 2.001188208394222
gamma_m 150.0 gamma_d 125.0 2.5372464543855386
gamma_m 150.0 gamma_d 150.0 nan
gamma_m 150.0 gamma_d 175.0 -2.5840841134673402
gamma_m 150.0 gamma_d 200.0 -2.0727783992021025
gamma_m 175.0 gamma_d 0.0 1.4645918875615234
gamma_m 175.0 gamma_d 25.0 1.4902570606397723
gamma_m 175.0 gamma_d 50.0 1.5585019217103644
gamma_m 175.0 gamma_d 75.0 1.6687875802778238
gamma_m 175.0 gamma_d 100.0 1.836572392913886
gamma_m 175.0 gamma_d 125.0 2.1098678647384412
gamma_m 175.0 gamma_d 150.0 2.6739764366927052
gamma_m 175.0 gamma_d 175.0 nan
gamma_m 175.0 gamma_d 200.0 -2.71637250521393
gamma_m 200.0 gamma_d 0.0 1.4645918875615231
gamma_m 200.0 gamma_d 25.0 1.4845441581729593
gamma_m 200.0 gamma_d 50.0 1.5377251238700236
gamma_m 200.0 gamma_d 75.0 1.6219368867750477
gamma_m 200.0 gamma_d 100.0 1.7436710272644398
gamma_m 200.0 gamma_d 125.0 1.9220789459321215
gamma_m 200.0 gamma_d 150.0 2.208757298511275
gamma_m 200.0 gamma_d 175.0 2.7980755038451366
gamma_m 200.0 gamma_d 200.0 nan
/tmp/ipykernel_172232/638675132.py:7: RuntimeWarning: invalid value encountered in double_scalars
  alpha = cube_root((np.pi*(gamma_m**2 - gamma_d**2)**3)/((2*gamma_d+gamma_m)**2 * (gamma_d-gamma_m)**4))
/tmp/ipykernel_172232/638675132.py:4: DeprecationWarning: In future, it will be an error for 'np.bool_' scalars to be interpreted as an index
  return math.pow(abs(i),float(1)/3) * (1,-1)[i<0]

3. implementation and analysis [36/42]

3.1. analytical derivations

3.1.1. thermodynamic spontaneity (global stability) [5/5]

  • I'm defining "thermodynamic spontaneity" as whether F(h=0) > F(h=1)
  • analyzed by plotting the phase boundary line across 2 parameters when F(h=0) == F(h=1)
  1. DONE setting up equality
    20210104_thermodynamic_spontaneity_kfix.png
    \begin{equation} \left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(3 V_{d}\right)^{2 / 3}+B=\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(3 V_d+\frac{S_i^{3/2}}{2 \sqrt{\pi }}\right)^{2/3} + \left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi+ P \eta-8 K C_{i} \sqrt{\pi S_{i}}+B \end{equation}
  2. DONE solving for \(K\)
    20210104_thermodynamic_k.png
    \begin{equation} \label{therm_spont_k} K = \frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}- P \eta}{8 \pi-8 C_{i} \sqrt{\pi S_{i}}} \end{equation}
  3. DONE solving for \(P\)
    20210104_thermodynamic_p.png
    \begin{equation} \label{therm_spont_p} P=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}-8 K \left(\pi - C_{i} \sqrt{\pi S_{i}}\right)}{\eta} \end{equation}
  4. DONE solving for \(\sigma\)
    20210104_thermodynamic_sigma.png
    \begin{equation} \label{therm_spont_sigma} \sigma=\frac{\frac{-\gamma_{m} S_{i}+8 K \left(\pi - C_{i} \sqrt{\pi S_{i}}\right) + P \eta}{\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}}-\gamma_{d} \epsilon+\gamma_{m} \beta}{\alpha+\beta-\frac{S_{i}}{\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}}} \end{equation}
  5. DONE solving for \(C_{i}\)
    20210104_thermodynamic_ci.png
    \begin{equation} \label{therm_spont_ci} C_{i} = \frac{\pi-\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \in\right)\left(\left(3 V_{d}\right)^{2 / 3}-\left(\frac{S_{i}^{3 / 2}}{2 \sqrt{\pi}}+3 V_{d}\right)^{2 / 3}\right)-\left(\sigma-\gamma_{m}\right) S_{i}-P_{\eta}}{8 K}}{\sqrt{\pi S_{i}}} \end{equation}
  6. any others?

3.1.2. kinetic spontaneity (local stability) [11/11]

  • I'm defining "kinetic spontaneity" as whether F'(h=0) is negative or positive
  1. DONE solving derivative
    20210105_f_derivative.png
    \begin{equation} \label{f-prime} F^{\prime}(h)=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right)\left(1-2 h^{2}\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(\frac{\left(\frac{3 h}{2}-h^{3}\right) S_{i}^{3 / 2}}{\sqrt{\pi}}+3 V_{d}\right)^{1 / 3}}-6 P \eta h^{2}+2\left(\left(\sigma-\gamma_{m}\right) S_{i}+8 K \pi\right) h+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation}
  2. DONE solving for h=0
    \begin{equation} \label{f-prime-alpha-beta} F^{\prime}(0)=\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi} \left(3 V_{d}\right)^{1 / 3}}+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation} \begin{equation} \label{f-prime} F^{\prime}(0)=\frac{\left(\gamma_{m} \alpha + \gamma_{d} \epsilon \right) S_{i}^{3 / 2}}{\sqrt{\pi} \left(3 V_{d}\right)^{1 / 3}}+3 P \eta-8 K C_{i} \sqrt{\pi S_{i}} \end{equation}
  3. DONE sign of droplet prefactor at h=0
    20210217_prefactor_sign.png

    the sign of the term concerning the droplet must always be positive.

  4. DONE solving for \(K\)
    20210106_kinetic_k.png
    \begin{equation} \label{kinetic_spont_k} K=\frac{\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+3 P \eta}{8 C_{i} \sqrt{\pi S_{i}}} \end{equation}
  5. DONE solving for \(P\)
    20210106_kinetic_p.png
    \begin{equation} \label{kinetic_spont_p} P=\frac{\frac{-\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}}{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+8 K C_{i} \sqrt{\pi S_{i}}}{3 \eta} \end{equation}
  6. DONE solving for \(\sigma\)
    20210106_kinetic_sigma.png
    \begin{equation} \label{kinetic_spont_sigma} \sigma=\frac{\frac{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)-\gamma_{d} \epsilon+\gamma_{m} \beta}{(\alpha+\beta)} \end{equation}
  7. DONE solving for \(C_{i}\)
    20210106_derivs_kinetic_ci.png
    \begin{equation} \label{kinetic_spont_ci} C_{i} = \frac{\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon\right) S_{i}^{3 / 2}} {\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}+3 P \eta}{8 K \sqrt{\pi S_{i}}} \end{equation}
  8. CANCELLED solving for \(\gamma_{m}\)
    20210106_derivs_kinetic_gammam.png
    \begin{equation} \label{kinetic_spont_gammam} \gamma_{m}=\frac{\frac{-\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)+\gamma_{d} \epsilon+\sigma \alpha+\sigma \beta}{\beta} \end{equation}
  9. CANCELLED solving for \(\gamma_{d}\)
    20210106_derivs_kinetic_gammad.png
    \begin{equation} \label{kinetic_spont_gammad} \gamma_{d}=\frac{\frac{\sqrt{\pi}\left(3 V_{d}\right)^{1 / 3}}{S_{i}^{3 / 2}}\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right)-\sigma \alpha-\left(\sigma-\gamma_{m}\right) \beta}{\epsilon} \end{equation}
  10. CANCELLED solving for \(S_{i}\)
  11. DONE solving for \(V_{d}\)
    20210106_derivs_kinetic_vd.png
    \begin{equation} \label{kinetic_spont_vd} V_{d}=\frac{\left(\frac{\left(\sigma \alpha+\left(\sigma-\gamma_{m}\right) \beta+\gamma_{d} \epsilon \right) S_{i}^{3 / 2}}{\left(8 K C_{i} \sqrt{\pi S_{i}}-3 P \eta \right) \sqrt{\pi}}\right)^{3}}{3} \end{equation}

3.2. python code [20/26]

3.2.1. configuration [12/13]

  1. import libraries, settings
    %gui qt
    import napari
    from decimal import Decimal
    import math
    import numpy as np
    import pandas as pd
    import matplotlib.pyplot as plt  # plotting
    import mpl_toolkits.mplot3d.axes3d as axes3d
    %matplotlib inline
    import seaborn as sns
    import scipy
    #sns.set(style='ticks')
    #sns.set_palette(sns.cubehelix_palette(10, start=.1, rot=-.75, reverse=True))
    sns.set_palette('hls',8)
    #sns.set_palette(sns.diverging_palette(250, 30, l=65, center="dark"))
    # sns.set_palette("coolwarm", 6)
    # from cycler import cycler
    SMALL_SIZE = 12
    MEDIUM_SIZE = 16
    BIGGER_SIZE = 20
    
    
    plt.rc('axes', titlesize=MEDIUM_SIZE)     # fontsize of the axes title
    plt.rc('axes', labelsize=MEDIUM_SIZE)    # fontsize of the x and y labels
    plt.rc('xtick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('ytick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('legend', fontsize=SMALL_SIZE)    # legend fontsize
    plt.rc('figure', titlesize=BIGGER_SIZE)  # fontsize of the figure title
    # default_cycler = (cycler(color=['indigo', 'darkgreen', 'purple', 'green', 'magenta', 'lime']))
    # plt.rc('axes', prop_cycle=default_cycler)
    
    pi = np.pi 
    
    #if kernel was started in notebook directory
    #fig_dir = 'figures/'
    #assuming kernel was started in home directory
    fig_dir = 'google_drive/grad_school/db_lab/code/analysis/llps_endocytosis_modeling/figures/'
    

    :results: nil:END:

    import IPython
    from tabulate import tabulate
    
    class OrgFormatter(IPython.core.formatters.BaseFormatter):
        def __call__(self, obj):
            try:
                return tabulate(obj, headers='keys',
                                tablefmt='orgtbl', showindex='always')
            except:
                return None
    
    ip = get_ipython()
    ip.display_formatter.formatters['text/org'] = OrgFormatter()
    

    :results: nil:END:

  2. set default fixed parameter values
    1. assigned here
      vesicle_radius = 50. #nm
      vesicle_diameter = 2*vesicle_radius #nm
      drop_radius = 1500. #nm
      S_i = 4.*pi*vesicle_radius**2 #nm^2
      #H_i_range = np.linspace(0.001*vesicle_radius,2*vesicle_radius,1000) #nm
      h_range = np.linspace(0.001,1,1000) #dimensionless 0-1
      H_i_range = h_range*(S_i/pi)**(1/2)
      sigma = 150. #arbitrary, N/nm^2
      A_t = 1000 #nm, arbitrary
      K = 100. #arbitrary, N/nm
      C_i = 1/vesicle_radius #nm^-1
      gamma_d = 200. #arbitrary N/nm^2
      gamma_m = 175. #arbitrary, N/nm^2
      gamma_i = gamma_m
      V_d = 4.*pi*drop_radius**2 #nm^3
      P = 1000/(10**9) #N/nm^3
      
      params = {'vesicle_radius' : vesicle_radius,
          'vesicle_diameter' : vesicle_diameter,
          'drop_radius' : drop_radius,
          'S_i' : S_i,
          'h_range' : h_range,
          'H_i_range' : H_i_range,
          'sigma' : sigma,
          'A_t' : A_t,
          'K' : K,
          'C_i' : C_i,
          'gamma_d' : gamma_d,
          'gamma_m' : gamma_m,
          'gamma_i' : gamma_i,
          'V_d' : V_d,
          'P':P}
      

      :results: nil:END:

    2. assigned from org table

      established in variables

      variable_df = pd.DataFrame(variable_list[1:], columns=variable_list[0])
      variable_df = variable_df.set_index('python variable')
      param_values = variable_df['value']
      variable_df
      
      python variable description value units reference type associated equation LaTeX variable
      Hi invagination height   \(nm\)   input (variable)   \(H_{i}\)
      sigma membrane tension 0.001 \(pN/nm\) hassinger_design_2017 fixed   \(\sigma\)
      At total plane radius 2000 \(nm\)   fixed   \(A_{t}\)
      vesicleradius final vesicle radius 50 \(nm\) kaksonen_mechanisms_2018 fixed   \(R_{v}\)
      K coat bending rigidity 160 \(pN nm\) dmitrieff_membrane_2015 fixed   \(K\)
      Ci spontaneous curvature 1/50 \(nm^{-1}\) hassinger_design_2017 fixed   \(C_{i}\)
      gammad droplet-cytosol surface tension 0.005 \(pN/nm\) jawerth_salt-dependent_2018 fixed   \(\gamma_{d}\)
      gammam droplet-membrane wetting tension 0.002 \(pN/nm\)   fixed   \(\gamma_{m}\)
        initial droplet plane radius       fixed    
      gammai droplet-invagination surface tension 0.002 \(pN/nm\)   fixed   \(\gamma_{i}\)
      P turgor pressure 0.1 \(pN/nm^2\) ma_mechanics_2019 fixed   \(P\)
      h relative invagination height   unitless   input (variable) rel-h \(h\)
      Vd droplet volume   \(nm^3\)   intermediate output   \(V_{d}\)
      Si invagination surface area   \(nm^2\)   intermediate output invagination-area \(S_{i}\)
      dropradius initial droplet spherical radius 600 \(nm\)   intermediate output    
      Sm droplet-flat membrane contact area   \(nm^2\)   intermediate output md-area \(S_{m}\)
      Ri invagination spherical radius   \(nm\)   intermediate output invagination-radius \(R_{i}\)
      Sd droplet-cytosol contact area   \(nm^2\)   intermediate output dc-contact-area \(S_{d}\)
      thetad droplet contact angle   rad   intermediate output droplet-contact-angle \(\theta_{d}\)
      thetai invagination contact angle   rad   intermediate output invagination-angle \(\theta_{i}\)
      Ai invagination plane radius   \(nm\)   intermediate output invagination-plane-radius \(A_{i}\)
      Ad droplet + invagination plane radius   \(nm\)   intermediate output droplet-plane-radius \(A_{d}\)
      So outside membrane surface area   \(nm^2\)   intermediate output outside-membrane-area \(S_{o}\)
      Vi invagination volume   \(nm^3\)   intermediate output invagination-volume \(V_{i}\)
      Vb droplet + invagination volume   \(nm^3\)   intermediate output di-volume \(V_{b}\)
      Hd droplet + invagination height   \(nm\)   intermediate output di-height \(H_{d}\)
      Rd droplet + invagination spherical radius   \(nm\)   intermediate output di-cap-radius \(R_{d}\)
      F free energy   \(pN \cdot nm\)   output f-forces f-shapes \(F\)
      membranetension membrane tension energy   \(pN \cdot nm\)   output    
      bendingenergy bending energy   \(pN \cdot nm\)   output    
      surfacetension surface tension energy   \(pN \cdot nm\)   output    
      wettingenergy wetting energy   \(pN \cdot nm\)   output    
      turgorpressure turgor pressure energy   \(pN \cdot nm\)   output    
      2÷Ri membrane curvature   \(nm^{-1}\)   output    
      vesicle_radius = param_values['vesicle_radius']
      vesicle_diameter = 2*vesicle_radius #nm
      drop_radius = param_values['drop_radius'] #nm
      S_i = 4.*pi*vesicle_radius**2 #nm^2
      #H_i_range = np.linspace(0.001*vesicle_radius,2*vesicle_radius,1000) #nm
      h_range = np.linspace(0.001,1,1000) #dimensionless 0-1
      H_i_range = h_range*(S_i/pi)**(1/2)
      sigma = param_values['sigma'] #arbitrary, N/nm^2
      A_t = param_values['A_t'] #nm, arbitrary
      K = param_values['K'] #arbitrary, N/nm
      C_i = 1/vesicle_radius #nm^-1
      gamma_d = param_values['gamma_d'] #arbitrary N/nm^2
      gamma_m = param_values['gamma_m'] #arbitrary, N/nm^2
      gamma_i = gamma_m
      V_d = 4.*pi*drop_radius**2 #nm^3
      P = param_values['P'] #N/nm^2
      
      params = {'vesicle_radius' : vesicle_radius,
          'vesicle_diameter' : vesicle_diameter,
          'drop_radius' : drop_radius,
          'S_i' : S_i,
          'h_range' : h_range,
          'H_i_range' : H_i_range,
          'sigma' : sigma,
          'A_t' : A_t,
          'K' : K,
          'C_i' : C_i,
          'gamma_d' : gamma_d,
          'gamma_m' : gamma_m,
          'gamma_i' : gamma_i,
          'V_d' : V_d,
          'P':P}
      

      :results: nil:END:

  3. define composite variables

    This is just to check what they are, I re-define them whenever I solve for F

    alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
    gamma_ratio = gamma_m/gamma_d
    alpha_refactored = np.cbrt((pi*(gamma_ratio+1)**3)/((gamma_ratio-1)*(gamma_ratio+2)**2))
    beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
    #epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
    epsilon = np.cbrt((8*pi)/((2+gamma_ratio)**2*(1-gamma_ratio)))
    eta = (S_i**(3/2))/(6*pi**(1/2))
    C = sigma*pi*A_t**2+S_i*(gamma_m-gamma_i)
    
    alpha, alpha_refactored, beta, epsilon, eta, C
    
    -1.3561885767612312 -1.356188576761231 1.3561885767612312 1.937412252516044 523598.7755982989 12566.370614359173
  4. custom functions [12/13]
    1. model functions
      1. pre-analytical version
        • [ ] rewrite each function variable as dictionary key (outpus['theta_d'] instead of theta_d for the first one) and just return outputs at end
          • should eliminate some redundancy

             def llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                           C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d,
                           P=P):
                 #H_i = h*(S_i/pi)**(1/2)
                 # might need to be -gamma_m/gamma_d
                 theta_d = np.arccos(gamma_m/gamma_d) #equation 8
                 #to make a full array out of the constant theta_d value
                 #theta_d = np.full_like(H_i, theta_d)  
                 A_i = ((S_i/pi)-H_i**2)**(1/2) #equation 7
                 V_i = (1/6)*pi*H_i*(3*A_i**2+H_i**2) #equation 6
             # this is actually wrong!
             #    R_i = (A_i**2-H_i**2)/(2*H_i) #equation 10
                 R_i = S_i/(2*pi*H_i) #equation 10 re-do
             # this is wrong!
             #    theta_i = np.arcsin(A_i/R_i) #equation 12
                 theta_i = np.arccos(1-S_i/(2*pi*R_i**2)) #equation 12 re-do
                 V_b = V_d+V_i # equation 5
            #      R_d_list = []
            #      for V in V_b:
            #          R_d_list.append(cube_root((3*V)/(pi*(2+np.cos(theta_d))*(1-np.cos(theta_d))**2))) 
                 R_d = ((3*V_b)/(pi*(2+np.cos(theta_d))*(1-np.cos(theta_d))**2))**(1/3) #equation 4
                 # R_d = np.array(R_d_list)
                 S_d = 2*pi*R_d**2*(1-np.cos(theta_d)) #equation 11
                 A_d = R_d*np.sin(theta_d) #equation 3
                 S_m = pi*(A_d**2-A_i**2) #equation 9
                 S_o = pi*(A_t**2-A_d**2) #equation 2
                 membrane_tension = sigma*(S_o+S_m+S_i)
                 bending_energy = (K/2)*S_i*((2/R_i)-(2*C_i))**2
                 surface_tension = gamma_d*S_d
                 wetting_energy = gamma_m*S_m+gamma_i*S_i
                 droplet_energy = surface_tension - wetting_energy
                 turgor_pressure = P*V_i
                 F = membrane_tension+bending_energy+surface_tension-wetting_energy+turgor_pressure #equation 1
            
                 outputs = {'theta_d':theta_d,'A_i':A_i,'V_i':V_i, 'R_i':R_i,
                            '2÷R_i':2/R_i,'theta_i':theta_i,'V_b':V_b,'R_d':R_d,
                           'S_d':S_d,'A_d':A_d,'S_m':S_m,'S_o':S_o,
                            'F':F,'membrane_tension':membrane_tension,
                            'bending_energy':bending_energy, 'bending_energy_flipped':-bending_energy,
                            'surface_tension':surface_tension,
                            'wetting_energy':-wetting_energy, 'wetting_energy_flipped':wetting_energy,
                            'droplet_energy':droplet_energy, 'turgor_pressure':turgor_pressure}
                 return outputs
            

            nil:END:

            plt.plot(H_i_range/100,llps_model()['F'])
            plt.xlabel('h')
            plt.ylabel('F')
            plt.tight_layout()
            plt.savefig(fig_dir+'llps_model_test.png')
            

            :results: c8fe3356fecff67fa663cdfd814c256ec1d29d36.png

      2. post-analytical version
        def llps_model_recbrt(h=h_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i,
                       V_d=V_d, P=P):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
            C = sigma*pi*A_t**2+S_i*(gamma_m-gamma_i)
        
            F = ((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)* np.cbrt(((S_i**(3/2))/(pi**(1/2)))*(3*h/2-h**3)+3*V_d)**2+ (sigma-gamma_m)*S_i*h**2+ 2*K*(C_i*S_i**(1/2)-2*pi**(1/2)*h)**2+ P*eta*(3*h-2*h**3)+ C )
        
        #    print(alpha, beta, epsilon, eta, C)
        
            return F
        

        nil:END:

        plt.plot(h_range,llps_model_recbrt())
        plt.xlabel('h')
        plt.ylabel('F')
        plt.tight_layout()
        plt.savefig(fig_dir+'llps_model_recbrt_test.png')
        

        :results: c8fe3356fecff67fa663cdfd814c256ec1d29d36.png

    2. general plotting function
      • figure out some way to get the right units for each thing being plotted
      def plot_simulation(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=None, variable_df=variable_df):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          y = outputs[output]
          plt.plot(x,y,label=label)
          plt.ylabel(variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
      

      :results: nil:END:

      1. testing
        plot_simulation()
        
        73a700d00e7b8c78670fcacedb4392900dd62848.png
    3. normalized plots
      def plot_simulation_normalized(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=None):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d)
          y_sub = outputs[output]-np.min(outputs[output])
          y_norm = y_sub/np.max(y_sub)
          plt.plot(x,y_norm,label=label)
          plt.ylabel(output)
          plt.xlabel('invagination distance/vesicle diameter')
      

      :results: nil:END:

    4. delta plots
      def plot_simulation_delta(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=None, variable_df=variable_df, color=None):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          output_0 = llps_model(H_i=H_i_range[0], sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
          y = outputs[output]-output_0[output]
          plt.plot(x,y,label=label, color=color)
          #plt.ylabel('∆ '+variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
          plt.xlabel('h: invagination distance/vesicle diameter')
      

      :results: nil:END:

    5. show geometry for given invagination distance
      def show_geometry(H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                        C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                        fig_width=10, fig_height=5,
                        num_plots=1, top=20, bottom=-250, left=-500, right=500):
          outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                   C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
      
          A_d = outputs['A_d']
          A_i = outputs['A_i']
          theta_d = outputs['theta_d']
          theta_i = outputs['theta_i']
          R_d = outputs['R_d']
          R_i = outputs['R_i']
      
          #fig, ax = plt.subplots(1,num_plots)
          #plt.subplot(1, num_plots, 1)
      
          fig = plt.figure()
          fig.set_size_inches((fig_width,fig_height))
          ax = fig.add_subplot(1, num_plots, 1)
      
          plt.hlines([0,0],[-A_t,A_d],[-A_d,A_t], color='r', label='$S_o$')
          plt.hlines([0,0],[-A_d,A_i],[-A_i,A_d], color='g', label='$S_m$')
      
          if V_d > 0:
            points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
            x_d = R_d*np.cos(points_d)
            y_d = R_d*np.sin(points_d)
            y_d_adj = y_d-np.min(y_d)
            plt.plot(x_d,-y_d_adj, color='b', label='$S_d$')
      
          points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
          x_i = R_i*np.cos(points_i)
          y_i = R_i*np.sin(points_i)
          y_i_adj = y_i-np.min(y_i)
          plt.plot(x_i,-y_i_adj, color='purple', label='$S_i$')
      
          ax.set_ylim(top=top, bottom=bottom)
          ax.set_xlim(left=left, right=right)
          ax.set_aspect('equal')
          ax.set(xlabel = 'distance along membrane (nm)', 
                 ylabel = 'distance into cytoplasm (nm)')
          plt.legend()
          plt.tight_layout()
      
          if num_plots > 1:
            return fig
      

      :results: nil:END:

      1. testing
        show_geometry(H_i=100)
        
        2d3471c918354444e5dbfed3090b9a23fc8c084c.png
    6. 3D geometry
      1. points
        %gui qt
        viewer = napari.Viewer(ndisplay=3)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_d, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_d
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='droplet', size=20, face_color='cyan',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_i, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_i
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='coat', size=20, face_color='purple',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        spherecoords = []
        for H_i in np.linspace(1,100,100):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            theta_m = 0.00000001
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 200), np.linspace(pi-theta_m, np.pi, 200)
            THETA, PHI = np.meshgrid(theta, phi)
            R = 0.1*A_t/np.sin(theta_m)
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)+0
            height = np.full_like(Z, H_i)
            spherecoords.append(np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T)
        
        spherecoords = np.concatenate(spherecoords)
        viewer.add_points(np.array(spherecoords), name='cell wall', size=20, face_color='red',
                          edge_color = 'transparent', blending='additive', opacity=0.6)
        
        <Points layer 'cell wall' at 0x7f3d216805b0>
        
        • I don't understand why this doesn't really work
        vertices = spherecoords
        #faces = np.array([[0, 0, 1, 2], [0, 1, 2, 3], [0, 2, 3, 4]])
        faces_list = []
        counter = 0
        for H_i in range(len(H_i_range)):
            for point in range(38):
                pointa = counter + point
                faces_list.append([pointa, pointa+1, pointa+2])
            counter += 40
        faces = np.array(faces_list)
        values = np.linspace(0, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer.add_surface(surface, colormap='viridis')  # add the surface
        
        <Surface layer 'surface' at 0x7f9405b20130>
        
      2. surface

        examples:

        vertices = np.array([[0, 0, 0, 0], [0, 10, 0, 20], [0, 20, 10, 0], [0, 0, 10, 10],
                             [1, 20, 0, 0], [1, 10, 0, 20], [1, 20, 10, 0], [1, 0, 10, 10]])
        faces = np.array([[0, 1, 2], [1, 2, 3],
                          [4,5,6],[5,6,7]])
        values = np.linspace(0.9, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer = napari.view_surface(surface)  # add the surface
        
        facesloop = []
        
        for t in (0,1):
            vert = vertices[vertices[:,0] == t]
            facest = scipy.spatial.Delaunay(vert[:,2:]).simplices + len(vert)*t
            facesloop.append(facest)
        
        facesd = np.concatenate(facesloop)
        
        viewer = napari.view_surface((vertices, facesd, values))
        
        coords = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T[:,1:]
        faces = scipy.spatial.Delaunay( coords[:,1:]).simplices
        values = np.linspace(0,1,len(coords))
        viewer = napari.view_surface((coords[:,1:], faces, values))
        
        %gui qt
        viewer = napari.Viewer(ndisplay=3)
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_d, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_d
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            #face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='droplet', colormap='green',
                           blending='translucent', shading='smooth')
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_i, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = R_i
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='coat', colormap='bop blue',
                           blending='additive', shading='smooth')
        
        spherecoords = []
        faces = []
        for H_i in range(1,101):
            outputs = llps_model(H_i=H_i, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
            A_d = outputs['A_d']
            A_i = outputs['A_i']
            theta_d = outputs['theta_d']
            theta_i = outputs['theta_i']
            theta_m = 0.00000001
            R_d = outputs['R_d']
            R_i = outputs['R_i']
            theta, phi = np.linspace(0,2* np.pi, 50), np.linspace(pi-theta_m, np.pi, 50)
            THETA, PHI = np.meshgrid(theta, phi)
            R = 0.1*A_t/np.sin(theta_m)
            X = R * np.sin(PHI) * np.cos(THETA)
            Y = R * np.sin(PHI) * np.sin(THETA)
            Z = R * np.cos(PHI)
            Z = Z-np.max(Z)+0
            height = np.full_like(Z, H_i)
            vertices = np.array([np.ravel(height), np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
            spherecoords.append(vertices)
        
            #face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            if H_i == 1:
                face = scipy.spatial.Delaunay(vertices[:,2:]).simplices + len(vertices)*(H_i-1)
            else:
                face = faces[-1] + len(vertices)
            faces.append(face)
        
        spherecoords = np.concatenate(spherecoords)
        faces = np.concatenate(faces)
        surface_values = np.full(len(spherecoords),1)
        viewer.add_surface((spherecoords,faces,surface_values), name='cell wall', colormap='bop blue',
                           blending='opaque', shading='smooth')
        
        def update_frame():
            """Update frame or time"""
            viewer.text_overlay.text = 'h = %1.2f\nblue: membrane\ngreen: condensate' % ((viewer.dims.current_step[0]+1)/100)
        
        viewer.dims.events.current_step.connect(lambda event: update_frame())
        
        viewer.text_overlay.visible=True
        viewer.text_overlay.font_size=20
        viewer.text_overlay.color = 'white'
        viewer.scale_bar.unit = 'nm'
        viewer.scale_bar.visible=True
        viewer.scale_bar.font_size=10
        

        nil:END:

        animation script

        from napari_animation import Animation
        
        animation = Animation(viewer)
        viewer.reset_view()
        viewer.dims.current_step = (0, 0,0,0)
        animation.capture_keyframe(steps=15)
        animation.capture_keyframe(steps=15)
        viewer.camera.zoom = 2.6
        viewer.camera.angles=(-0.026467524771271397, 0.5719902891444765, 11.550928988495166)
        animation.capture_keyframe(steps=50)
        viewer.dims.current_step = (viewer.dims.range[0][1]-1, 0,0,0)
        animation.capture_keyframe(steps=int(viewer.dims.range[0][1]-1))
        viewer.reset_view()
        animation.capture_keyframe(steps=50)
        viewer.dims.current_step = (0, 0,0,0)
        animation.capture_keyframe(steps=int(viewer.dims.range[0][1]-1))
        animation.animate('~/google_drive/grad_school/db_lab/code/analysis/llps_endocytosis_modeling/scripted_surface_animation.mov')
        

        :results:

        IMAGEIO FFMPEG_WRITER WARNING: input image is not divisible by macro_block_size=16, resizing from (2292, 1318) to (2304, 1328) to ensure video compatibility with most codecs and players. To prevent resizing, make your input image divisible by the macro_block_size or set the macro_block_size to 1 (risking incompatibility).
        Rendering frame  1 of 316
        Rendering frame  2 of 316
        Rendering frame  3 of 316
        Rendering frame  4 of 316
        Rendering frame  5 of 316
        Rendering frame  6 of 316
        Rendering frame  7 of 316
        Rendering frame  8 of 316
        Rendering frame  9 of 316
        Rendering frame  10 of 316
        Rendering frame  11 of 316
        Rendering frame  12 of 316
        Rendering frame  13 of 316
        Rendering frame  14 of 316
        Rendering frame  15 of 316
        Rendering frame  16 of 316
        Rendering frame  17 of 316
        Rendering frame  18 of 316
        Rendering frame  19 of 316
        Rendering frame  20 of 316
        Rendering frame  21 of 316
        Rendering frame  22 of 316
        Rendering frame  23 of 316
        Rendering frame  24 of 316
        Rendering frame  25 of 316
        Rendering frame  26 of 316
        Rendering frame  27 of 316
        Rendering frame  28 of 316
        Rendering frame  29 of 316
        Rendering frame  30 of 316
        Rendering frame  31 of 316
        Rendering frame  32 of 316
        Rendering frame  33 of 316
        Rendering frame  34 of 316
        Rendering frame  35 of 316
        Rendering frame  36 of 316
        Rendering frame  37 of 316
        Rendering frame  38 of 316
        Rendering frame  39 of 316
        Rendering frame  40 of 316
        Rendering frame  41 of 316
        Rendering frame  42 of 316
        Rendering frame  43 of 316
        Rendering frame  44 of 316
        Rendering frame  45 of 316
        Rendering frame  46 of 316
        Rendering frame  47 of 316
        Rendering frame  48 of 316
        Rendering frame  49 of 316
        Rendering frame  50 of 316
        Rendering frame  51 of 316
        Rendering frame  52 of 316
        Rendering frame  53 of 316
        Rendering frame  54 of 316
        Rendering frame  55 of 316
        Rendering frame  56 of 316
        Rendering frame  57 of 316
        Rendering frame  58 of 316
        Rendering frame  59 of 316
        Rendering frame  60 of 316
        Rendering frame  61 of 316
        Rendering frame  62 of 316
        Rendering frame  63 of 316
        Rendering frame  64 of 316
        Rendering frame  65 of 316
        Rendering frame  66 of 316
        Rendering frame  67 of 316
        Rendering frame  68 of 316
        Rendering frame  69 of 316
        Rendering frame  70 of 316
        Rendering frame  71 of 316
        Rendering frame  72 of 316
        Rendering frame  73 of 316
        Rendering frame  74 of 316
        Rendering frame  75 of 316
        Rendering frame  76 of 316
        Rendering frame  77 of 316
        Rendering frame  78 of 316
        Rendering frame  79 of 316
        Rendering frame  80 of 316
        Rendering frame  81 of 316
        Rendering frame  82 of 316
        Rendering frame  83 of 316
        Rendering frame  84 of 316
        Rendering frame  85 of 316
        Rendering frame  86 of 316
        Rendering frame  87 of 316
        Rendering frame  88 of 316
        Rendering frame  89 of 316
        Rendering frame  90 of 316
        Rendering frame  91 of 316
        Rendering frame  92 of 316
        Rendering frame  93 of 316
        Rendering frame  94 of 316
        Rendering frame  95 of 316
        Rendering frame  96 of 316
        Rendering frame  97 of 316
        Rendering frame  98 of 316
        Rendering frame  99 of 316
        Rendering frame  100 of 316
        Rendering frame  101 of 316
        Rendering frame  102 of 316
        Rendering frame  103 of 316
        Rendering frame  104 of 316
        Rendering frame  105 of 316
        Rendering frame  106 of 316
        Rendering frame  107 of 316
        Rendering frame  108 of 316
        Rendering frame  109 of 316
        Rendering frame  110 of 316
        Rendering frame  111 of 316
        Rendering frame  112 of 316
        Rendering frame  113 of 316
        Rendering frame  114 of 316
        Rendering frame  115 of 316
        Rendering frame  116 of 316
        Rendering frame  117 of 316
        Rendering frame  118 of 316
        Rendering frame  119 of 316
        Rendering frame  120 of 316
        Rendering frame  121 of 316
        Rendering frame  122 of 316
        Rendering frame  123 of 316
        Rendering frame  124 of 316
        Rendering frame  125 of 316
        Rendering frame  126 of 316
        Rendering frame  127 of 316
        Rendering frame  128 of 316
        Rendering frame  129 of 316
        Rendering frame  130 of 316
        Rendering frame  131 of 316
        Rendering frame  132 of 316
        Rendering frame  133 of 316
        Rendering frame  134 of 316
        Rendering frame  135 of 316
        Rendering frame  136 of 316
        Rendering frame  137 of 316
        Rendering frame  138 of 316
        Rendering frame  139 of 316
        Rendering frame  140 of 316
        Rendering frame  141 of 316
        Rendering frame  142 of 316
        Rendering frame  143 of 316
        Rendering frame  144 of 316
        Rendering frame  145 of 316
        Rendering frame  146 of 316
        Rendering frame  147 of 316
        Rendering frame  148 of 316
        Rendering frame  149 of 316
        Rendering frame  150 of 316
        Rendering frame  151 of 316
        Rendering frame  152 of 316
        Rendering frame  153 of 316
        Rendering frame  154 of 316
        Rendering frame  155 of 316
        Rendering frame  156 of 316
        Rendering frame  157 of 316
        Rendering frame  158 of 316
        Rendering frame  159 of 316
        Rendering frame  160 of 316
        Rendering frame  161 of 316
        Rendering frame  162 of 316
        Rendering frame  163 of 316
        Rendering frame  164 of 316
        Rendering frame  165 of 316
        Rendering frame  166 of 316
        Rendering frame  167 of 316
        Rendering frame  168 of 316
        Rendering frame  169 of 316
        Rendering frame  170 of 316
        Rendering frame  171 of 316
        Rendering frame  172 of 316
        Rendering frame  173 of 316
        Rendering frame  174 of 316
        Rendering frame  175 of 316
        Rendering frame  176 of 316
        Rendering frame  177 of 316
        Rendering frame  178 of 316
        Rendering frame  179 of 316
        Rendering frame  180 of 316
        Rendering frame  181 of 316
        Rendering frame  182 of 316
        Rendering frame  183 of 316
        Rendering frame  184 of 316
        Rendering frame  185 of 316
        Rendering frame  186 of 316
        Rendering frame  187 of 316
        Rendering frame  188 of 316
        Rendering frame  189 of 316
        Rendering frame  190 of 316
        Rendering frame  191 of 316
        Rendering frame  192 of 316
        Rendering frame  193 of 316
        Rendering frame  194 of 316
        Rendering frame  195 of 316
        Rendering frame  196 of 316
        Rendering frame  197 of 316
        Rendering frame  198 of 316
        Rendering frame  199 of 316
        Rendering frame  200 of 316
        Rendering frame  201 of 316
        Rendering frame  202 of 316
        Rendering frame  203 of 316
        Rendering frame  204 of 316
        Rendering frame  205 of 316
        Rendering frame  206 of 316
        Rendering frame  207 of 316
        Rendering frame  208 of 316
        Rendering frame  209 of 316
        Rendering frame  210 of 316
        Rendering frame  211 of 316
        Rendering frame  212 of 316
        Rendering frame  213 of 316
        Rendering frame  214 of 316
        Rendering frame  215 of 316
        Rendering frame  216 of 316
        Rendering frame  217 of 316
        Rendering frame  218 of 316
        Rendering frame  219 of 316
        Rendering frame  220 of 316
        Rendering frame  221 of 316
        Rendering frame  222 of 316
        Rendering frame  223 of 316
        Rendering frame  224 of 316
        Rendering frame  225 of 316
        Rendering frame  226 of 316
        Rendering frame  227 of 316
        Rendering frame  228 of 316
        Rendering frame  229 of 316
        Rendering frame  230 of 316
        Rendering frame  231 of 316
        Rendering frame  232 of 316
        Rendering frame  233 of 316
        Rendering frame  234 of 316
        Rendering frame  235 of 316
        Rendering frame  236 of 316
        Rendering frame  237 of 316
        Rendering frame  238 of 316
        Rendering frame  239 of 316
        Rendering frame  240 of 316
        Rendering frame  241 of 316
        Rendering frame  242 of 316
        Rendering frame  243 of 316
        Rendering frame  244 of 316
        Rendering frame  245 of 316
        Rendering frame  246 of 316
        Rendering frame  247 of 316
        Rendering frame  248 of 316
        Rendering frame  249 of 316
        Rendering frame  250 of 316
        Rendering frame  251 of 316
        Rendering frame  252 of 316
        Rendering frame  253 of 316
        Rendering frame  254 of 316
        Rendering frame  255 of 316
        Rendering frame  256 of 316
        Rendering frame  257 of 316
        Rendering frame  258 of 316
        Rendering frame  259 of 316
        Rendering frame  260 of 316
        Rendering frame  261 of 316
        Rendering frame  262 of 316
        Rendering frame  263 of 316
        Rendering frame  264 of 316
        Rendering frame  265 of 316
        Rendering frame  266 of 316
        Rendering frame  267 of 316
        Rendering frame  268 of 316
        Rendering frame  269 of 316
        Rendering frame  270 of 316
        Rendering frame  271 of 316
        Rendering frame  272 of 316
        Rendering frame  273 of 316
        Rendering frame  274 of 316
        Rendering frame  275 of 316
        Rendering frame  276 of 316
        Rendering frame  277 of 316
        Rendering frame  278 of 316
        Rendering frame  279 of 316
        Rendering frame  280 of 316
        Rendering frame  281 of 316
        Rendering frame  282 of 316
        Rendering frame  283 of 316
        Rendering frame  284 of 316
        Rendering frame  285 of 316
        Rendering frame  286 of 316
        Rendering frame  287 of 316
        Rendering frame  288 of 316
        Rendering frame  289 of 316
        Rendering frame  290 of 316
        Rendering frame  291 of 316
        Rendering frame  292 of 316
        Rendering frame  293 of 316
        Rendering frame  294 of 316
        Rendering frame  295 of 316
        Rendering frame  296 of 316
        Rendering frame  297 of 316
        Rendering frame  298 of 316
        Rendering frame  299 of 316
        Rendering frame  300 of 316
        Rendering frame  301 of 316
        Rendering frame  302 of 316
        Rendering frame  303 of 316
        Rendering frame  304 of 316
        Rendering frame  305 of 316
        Rendering frame  306 of 316
        Rendering frame  307 of 316
        Rendering frame  308 of 316
        Rendering frame  309 of 316
        Rendering frame  310 of 316
        Rendering frame  311 of 316
        Rendering frame  312 of 316
        Rendering frame  313 of 316
        Rendering frame  314 of 316
        Rendering frame  315 of 316
        Rendering frame  316 of 316
        
      3. probably don't need this stuff
        num_plots = 1
        outputs = llps_model(H_i=90, sigma=sigma, A_t=A_t, S_i=S_i, K=K, C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P)
        
        A_d = outputs['A_d']
        A_i = outputs['A_i']
        theta_d = outputs['theta_d']
        theta_i = outputs['theta_i']
        R_d = outputs['R_d']
        R_i = outputs['R_i']
        
        #fig, ax = plt.subplots(1,num_plots)
        #plt.subplot(1, num_plots, 1)
        
        #fig = plt.figure()
        #fig.set_size_inches((fig_width,fig_height))
        #ax = fig.add_subplot(1, num_plots, 1)
        
        #plt.hlines([0,0],[-A_t,A_d],[-A_d,A_t], color='r', label='$S_o$')
        #plt.hlines([0,0],[-A_d,A_i],[-A_i,A_d], color='g', label='$S_m$')
        
        if V_d > 0:
            points_d = np.linspace(pi/2-theta_d, theta_d+pi/2, 1000)
            x_d = R_d*np.cos(points_d)
            y_d = R_d*np.sin(points_d)
            y_d_adj = y_d-np.min(y_d)
        
        points_i = np.linspace(pi/2-theta_i, theta_i+pi/2, 1000)
        x_i = R_i*np.cos(points_i)
        y_i = R_i*np.sin(points_i)
        y_i_adj = y_i-np.min(y_i)
        
        xy_i = np.array([y_i,x_i]).T
        viewer = napari.Viewer()
        viewer.add_points(xy_i)
        
        <Points layer 'xy_i' at 0x7f9405887370>
        
        nx, ny, nz = (3, 3, 4)
        x = np.linspace(0, 10, nx)
        y = np.linspace(0, 10, ny)
        z = np.linspace(0, 10, nz)
        
        xv, yv, zv = np.meshgrid(x, y, z)
        zyxv = np.array([np.ravel(zv), np.ravel(yv), np.ravel(xv)]).T
        zyxv
        
        viewer = napari.Viewer()
        viewer.add_points(zyxv)
        
        <Points layer 'zyxv' at 0x7f93f2b41fd0>
        
        import numpy as np
        import matplotlib.pyplot as plt
        from mpl_toolkits import mplot3d
        
        ax = plt.axes(projection="3d")
        xlist=np.linspace(-1.0,1.0,10)*100
        ylist=np.linspace(-1.0,1.0,10)*100
        r=np.linspace(1.0,1.0,10)*100
        X,Y= np.meshgrid(xlist,ylist)
        
        Z=np.sqrt(r**2-X**2-Y**2) #Use np.sqrt like you had before
        
        spherecoords = np.array([np.ravel(Z), np.ravel(Y), np.ravel(X)]).T
        viewer.add_points(spherecoords)
        
        vertices = spherecoords
        faces = np.array([[0, 50, 99], [1, 2, 3]])
        values = np.linspace(1, 1, len(vertices))
        surface = (vertices, faces, values)
        
        viewer.add_surface(surface)
        
        cp=ax.plot_wireframe(X,Y,Z,color="r")
        cp=ax.plot_wireframe(X,Y,-Z,color="r") # Now plot the bottom half
        
        
        plt.title('2D Contour Plot of The unit sphere')
        plt.show()
        
        <ipython-input-21-c34f0052e1a6>:11: RuntimeWarning: invalid value encountered in sqrt
          Z=np.sqrt(r**2-X**2-Y**2) #Use np.sqrt like you had before
        
        e0cc5e561ea03e27b011a913d05aff6ff077c5fd.png
    7. side-by-side delta plot and equilibrium geometry
      def plot_simulation_geometry(x=H_i_range/vesicle_diameter, output='F',
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       fig_width=15, fig_height=5, label=None, variable_df=variable_df,
                      num_plots=2, top=20, bottom=-250, left=-500, right=500):
          outputs = llps_model(H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                    C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d)
          min_energy_ix = np.argmin(outputs['F'])
          H_i_eq = H_i_range[min_energy_ix]
          fig = show_geometry(H_i=H_i_eq, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                         C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                         fig_width=fig_width, fig_height=fig_height,
                         num_plots=num_plots, top=top, bottom=bottom, left=left, right=right)
          plt.title('geometry at equilibrium')
          ax2 = fig.add_subplot(1,2,2)
          #y = outputs[output]
          #plt.plot(x,y,label=label)
          #plt.ylabel(output)
      
          plot_simulation_delta(x=x, output=output,
                       H_i=H_i_range, sigma=sigma, A_t=A_t, S_i=S_i, K=K,
                       C_i=C_i, gamma_d=gamma_d, gamma_m=gamma_m, gamma_i=gamma_i, V_d=V_d, P=P,
                       label=label, variable_df=variable_df)
      
          plt.tight_layout()
      

      :results: nil:END:

      1. testing
        plot_simulation_geometry()
        
        c061a48c607c1e529382b8ea7629e0dacb2719aa.png
    8. side-by-side F plots and parameter phase diagram
      def boundary_test(x_param, y_param, x_range, boundary, params=params,
                        factors=[0.7,1,1.3], y_offset=0.2, x_scale='log'):
          fparams = params.copy()
          fig, (ax1, ax2)  = plt.subplots(nrows=1, ncols=2, figsize=(10,4))
      
          if x_scale == 'log':
              ax1.semilogx(x_range, boundary, color='black')
          else:
              ax1.plot(x_range, boundary, color='black')
      
          for factor in factors:
              middle = len(x_range)//2
              x_pick = x_range[middle]*factor
              boundary_pick = boundary[middle]+(factor-1)*np.sqrt(boundary[middle]**2)*y_offset
              fparams[x_param] = x_pick
              fparams[y_param] = boundary_pick
              F = llps_model_recbrt(h=fparams['h_range'], sigma=fparams['sigma'], A_t=fparams['A_t'],
                      S_i=fparams['S_i'], K=fparams['K'], C_i=fparams['C_i'], gamma_d=fparams['gamma_d'],
                      gamma_m=fparams['gamma_m'], gamma_i=fparams['gamma_i'], V_d=fparams['V_d'], P=fparams['P'])
              delta_F = F-F[0]
      
              ax1.scatter(x_pick, boundary_pick)
              ax2.plot(h_range,delta_F)
      
          ax2.hlines(0, 0,1.05, color='black')
      
          ax1.set(xlabel=x_param, ylabel=y_param)
          ax2.set(xlabel='h', ylabel='∆F')
          plt.tight_layout()
      

      :results: nil:END:

    9. global stability parameter analysis [4/4]
      1. DONE \(K\)

        from therm_spont_k

        def thermo_spont_K(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, P=P, C_i=C_i):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            K = ((-(sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*
                  (np.cbrt(3*V_d+(S_i**(3/2))/(2*pi**(1/2)))**2-
                    np.cbrt(3*V_d)**2)+2*P*eta-(sigma-gamma_m)*S_i-
                  3*P*eta)/(8*pi-8*C_i*(pi*S_i)**(1/2)))
        
            return K
        

        :results: nil:END:

        1. testing
          V_d_range = np.logspace(6,9)
          K_boundary = thermo_spont_K(V_d=V_d_range)
          boundary_test(x_param='V_d', y_param='K', x_range=V_d_range, boundary=K_boundary)
          
          f2003c3926174656e20a4b3bc9d05d90a3d8e8ed.png
          C_i_range = np.logspace(-4,-2.1,1000)
          K_boundary = thermo_spont_K(C_i=C_i_range)
          boundary_test(x_param='C_i', y_param='K', x_range=C_i_range, boundary=K_boundary,
                        y_offset=-0.3)
          
          1420f241c02fa916e57a933e084cca10a86bb665.png
          • this kind of makes sense - the lower the spontaneous curvature (flatter), the more the energetics are sensitive to rigidity
          x_range = np.logspace(6,9)
          thermo_boundary = thermo_spont_K(V_d=x_range)
          kinetic_boundary = kinetic_spont_K(V_d=x_range)
          both_boundary_test(x_param='V_d', y_param='K', x_range=x_range,
                             thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                             y_offset=2)
          
          ---------------------------------------------------------------------------
          NameError                                 Traceback (most recent call last)
          Cell In [24], line 3
                1 x_range = np.logspace(6,9)
                2 thermo_boundary = thermo_spont_K(V_d=x_range)
          ----> 3 kinetic_boundary = kinetic_spont_K(V_d=x_range)
                4 both_boundary_test(x_param='V_d', y_param='K', x_range=x_range,
                5                    thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                6                    y_offset=2)
          
          NameError: name 'kinetic_spont_K' is not defined
          
      2. DONE \(P\)
        def thermo_spont_P(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            P = ((sigma*alpha + (sigma - gamma_m)*beta + gamma_d*epsilon)*(np.cbrt(3*V_d)**2 -
                 np.cbrt(S_i**(3/2)/(2*np.sqrt(pi)) + 3*V_d)**2) -
                 (sigma-gamma_m)*S_i - 8*K*(pi - C_i*np.sqrt(pi*S_i)))/ eta
        
            return P
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(6,9)
          boundary = thermo_spont_P(V_d=x_range)
          boundary_test(x_param='V_d', y_param='P', x_range=x_range, boundary=boundary,
                        y_offset=-0.3)
          
          e1f50629882f81b82adf66ee1e8b21af1629cd8e.png
      3. DONE \(\sigma\)
        def thermo_spont_sigma(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            sigma = ((-gamma_m*S_i+8*K*(pi-C_i*np.sqrt(pi*S_i))+P*eta)/
                    ((3*V_d)**(2/3)-((S_i**(3/2))/
                        (2*np.sqrt(pi))+3*V_d)**(2/3))-
                     gamma_d*epsilon+gamma_m*beta)/(alpha+beta-
                    (S_i)/((3*V_d)**(2/3)-((S_i**(3/2))/
                            (2*np.sqrt(pi))+3*V_d)**(2/3)))
        
            return sigma
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(6,9)
          boundary = thermo_spont_sigma(V_d=x_range)
          boundary_test(x_param='V_d', y_param='sigma', x_range=x_range, boundary=boundary,
                        y_offset=-0.02)
          
          e7b526abe35c6ee6ee4e49a8a049fbfd1919a276.png
      4. DONE \(C_{i}\)
        def thermo_spont_C_i(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, sigma=sigma):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            C_i = (pi-((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*(np.cbrt(3*V_d)**2-np.cbrt((S_i**(3/2))/(2*np.sqrt(pi))+3*V_d)**2)-(sigma-gamma_m)*S_i-P*eta)/(8*K))/(np.sqrt(pi*S_i))
        
            return C_i
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(2,7)
          boundary = thermo_spont_C_i(V_d=x_range)
          boundary_test(x_param='V_d', y_param='C_i', x_range=x_range, boundary=boundary,
                        y_offset=0.2)
          
          26c718686e40c5395d5dbc51cf555d5602f0d796.png
          • it doesn't make sense that Ci would ever be negative…
    10. local stability parameter analysis [8/8]
      1. DONE \(K\)

        from kinetic_spont_k

        def kinetic_spont_K(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, P=P, C_i=C_i):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            K = (((sigma*alpha+(sigma-gamma_m)*beta+gamma_d*epsilon)*S_i**(3/2))/(np.sqrt(pi)*np.cbrt(3*V_d)) + 3*P*eta)/(8*C_i*np.sqrt(pi*S_i))
        
            return K
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_K(V_d=x_range)
          boundary_test(x_param='V_d', y_param='K', x_range=x_range, boundary=boundary,
                        y_offset=0.5)
          
          fa622209b5fd08a71c1fd5e6d57832eb4c129095.png
      2. DONE \(P\)

        from kinetic_spont_p

        def kinetic_spont_P(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, C_i=C_i):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            P = (-((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 (np.sqrt(pi)*np.cbrt(3*V_d)) + 8*K*C_i*np.sqrt(pi*S_i))/(3*eta)
        
            return P
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_P(V_d=x_range)
          boundary_test(x_param='V_d', y_param='P', x_range=x_range, boundary=boundary,
                        y_offset=-0.5)
          
          0af00632ca384881bfaef3c50556abefa0833290.png
      3. DONE \(\sigma\)

        nil:END: :results: nil:END: :results: nil:END:

        :RESULTS:

        ---------------------------------------------------------------------------
        NameError                                 Traceback (most recent call last)
        <ipython-input-35-594c148ecf9a> in <module>
        ----> 1 def kinetic_spont_sigma(P=P, gamma_m=gamma_m, gamma_d=gamma_d,
              2                   V_d=V_d, S_i=S_i, K=K, C_i=C_i):
              3 
              4     alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
              5     beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
        
        NameError: name 'P' is not defined
        
        1. testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_sigma(V_d=x_range)
          boundary_test(x_param='V_d', y_param='sigma', x_range=x_range, boundary=boundary,
                        y_offset=-0.5)
          
          ---------------------------------------------------------------------------
          TypeError                                 Traceback (most recent call last)
          Cell In [36], line 2
                1 x_range = np.logspace(6,9)
          ----> 2 boundary = kinetic_spont_sigma(V_d=x_range)
                3 boundary_test(x_param='V_d', y_param='sigma', x_range=x_range, boundary=boundary,
                4               y_offset=-0.5)
          
          Cell In [35], line 9, in kinetic_spont_sigma(P, gamma_m, gamma_d, V_d, S_i, K, C_i)
                6 epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
                7 eta = (S_i**(3/2))/(6*pi**(1/2))
          ----> 9 sigma = (( np.sqrt(pi)/S_i**(3/2) )*np.cbrt(3*V_d)(8*K*C_i*np.sqrt(pi*S_i)-3*P*eta)-
               10          gamma_d*epsilon+gamma_m*beta)/(alpha+beta)
               12 return sigma
          
          TypeError: 'numpy.ndarray' object is not callable
          

          :END:

      4. DONE \(C_{i}\)

        from kinetic_spont_ci

        def kinetic_spont_C_i(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          V_d=V_d, S_i=S_i, K=K, P=P):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            C_i = (((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 (np.sqrt(pi)*np.cbrt(3*V_d)) + 3*P*eta)/(8*K*np.sqrt(pi*S_i))
        
            return C_i
        

        :results: nil:END:

        1. testing
          x_range = np.logspace(6,9)
          boundary = kinetic_spont_C_i(V_d=x_range)
          boundary_test(x_param='V_d', y_param='C_i', x_range=x_range, boundary=boundary,
                        y_offset=0.5)
          
          d6c895cc30f3fee57f38868f32cc1510dac66fc3.png
      5. CANCELLED \(\gamma_m\)
      6. CANCELLED \(\gamma_d\)
      7. CANCELLED \(S_{i}\)
      8. DONE \(V_{d}\)

        from kinetic_spont_ci

        def kinetic_spont_V_d(sigma=sigma, gamma_m=gamma_m, gamma_d=gamma_d,
                          C_i=C_i, S_i=S_i, K=K, P=P):
        
            alpha = np.cbrt((pi*(gamma_m**2-gamma_d**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            beta = np.cbrt((pi*(gamma_d**2-gamma_m**2)**3)/((2*gamma_d+gamma_m)**2*(gamma_d-gamma_m)**4))
            epsilon = np.cbrt((8*pi)/((2+gamma_m/gamma_d)**2*(1-gamma_m/gamma_d)))
            eta = (S_i**(3/2))/(6*pi**(1/2))
        
            V_d = (1/3)*(((sigma*alpha + (sigma - gamma_m)*beta +gamma_d*epsilon)*S_i**(3/2))/
                 ((8*K*C_i*np.sqrt(pi*S_i)-3*P*eta)*np.sqrt(pi)))**3
        
            return V_d
        

        :results: nil:END:

        1. testing
          x_range = np.linspace(0.0002, 0.002)
          boundary = kinetic_spont_V_d(P=x_range)
          boundary_test(x_param='P', y_param='V_d', x_range=x_range, boundary=boundary,
                        y_offset=-5, x_scale='linear')
          
          47b70c0c626b8d536e2e83e12b8f0aed54dddab4.png
    11. side-by-side F plots with both phase boundaries
      def both_boundary_test(x_param, y_param, x_range, thermo_boundary, kinetic_boundary,
                             params=params, factors=[0.7,1,1.3], y_offset=0.2, x_scale='log', y_scale='linear',
                             variable_df=variable_df):
          fparams = params.copy()
          fig, (ax1, ax2)  = plt.subplots(nrows=1, ncols=2, figsize=(10,4))
      
      
          # if x_scale == 'log':
          #     ax1.semilogx(x_range, thermo_boundary, color='black', label='thermodynamic')
          #     ax1.semilogx(x_range, kinetic_boundary, color='blue', label='kinetic')
          # # else:
          #     ax1.plot(x_range, boundary, color='black')
      
          if y_offset > 0:
              fill_range = np.max(kinetic_boundary)
          if y_offset < 0:
              fill_range = np.min(kinetic_boundary)
      
          ax1.loglog(x_range, thermo_boundary, color='yellow')
          ax1.fill_between(x_range, y1=thermo_boundary, y2=fill_range*1.3, color='yellow', alpha=0.8, label='globally stable')
          ax1.loglog(x_range, kinetic_boundary, color='blue')
          ax1.fill_between(x_range, y1=kinetic_boundary, y2=fill_range*1.3, color='blue', alpha=0.8, label='locally stable')
      
          ax1.legend()
      
          for i, factor in enumerate(factors):
              middle = len(x_range)//2
              x_pick = x_range[middle]*factor
              for j, boundary in enumerate([thermo_boundary, kinetic_boundary]):
                boundary_pick = boundary[middle]+(factor-1)*np.sqrt(thermo_boundary[middle]**2)*y_offset
                fparams[x_param] = x_pick
                fparams[y_param] = boundary_pick
                F = llps_model_recbrt(h=fparams['h_range'], sigma=fparams['sigma'], A_t=fparams['A_t'],
                        S_i=fparams['S_i'], K=fparams['K'], C_i=fparams['C_i'], gamma_d=fparams['gamma_d'],
                        gamma_m=fparams['gamma_m'], gamma_i=fparams['gamma_i'], V_d=fparams['V_d'], P=fparams['P'])
                delta_F = F-F[0]
      
                # ax1.scatter(x_pick, boundary_pick, linewidths=5, color=((j)/2, (j)/2, (j+1)/2, (i+1)/len(factors)))
                # ax2.plot(h_range,delta_F, color=((j)/2, (j)/2, (j+1)/2, (i+1)/len(factors)))
      
                ax1.scatter(x_pick, boundary_pick, linewidths=5)
                ax2.plot(h_range,delta_F)
      
          ax2.hlines(0, -0.05,1.05, color='black')
      
          ax1.set(xlabel=variable_df.loc[x_param]['description']+' ('+variable_df.loc[x_param]['units']+')',
                  ylabel=variable_df.loc[y_param]['description']+' ('+variable_df.loc[y_param]['units']+')')
      #            xscale=x_scale)
          ax1.set_xscale(x_scale)
          ax1.set_yscale(y_scale)
          ax2.set(xlabel='h: relative invagination height', ylabel='∆F ($pN \cdot nm$)')
          plt.tight_layout()
      

      :results: nil:END:

3.2.2. result plots [8/9]

  1. plot all outputs vs. invagination height for default parameters
    1. as individual plots
      outputs = llps_model()
      outputs.pop('theta_d')
      # outputs.pop('2÷R_i')
      variables = list(outputs.keys())
      
      for variable in variables:
          plt.figure()
          plot_simulation(output=variable)
          plt.xlabel('h: invagination distance/vesicle diameter')
          plt.savefig(fig_dir+'%s-vs-H_i.png' %(variable))
      
      ---------------------------------------------------------------------------
      KeyError                                  Traceback (most recent call last)
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexes/base.py:3803, in Index.get_loc(self, key, method, tolerance)
         3802 try:
      -> 3803     return self._engine.get_loc(casted_key)
         3804 except KeyError as err:
      
      File ~/.local/lib/python3.10/site-packages/pandas/_libs/index.pyx:138, in pandas._libs.index.IndexEngine.get_loc()
      
      File ~/.local/lib/python3.10/site-packages/pandas/_libs/index.pyx:165, in pandas._libs.index.IndexEngine.get_loc()
      
      File pandas/_libs/hashtable_class_helper.pxi:5745, in pandas._libs.hashtable.PyObjectHashTable.get_item()
      
      File pandas/_libs/hashtable_class_helper.pxi:5753, in pandas._libs.hashtable.PyObjectHashTable.get_item()
      
      KeyError: 'bending_energy_flipped'
      
      The above exception was the direct cause of the following exception:
      
      KeyError                                  Traceback (most recent call last)
      Cell In [42], line 8
            6 for variable in variables:
            7     plt.figure()
      ----> 8     plot_simulation(output=variable)
            9     plt.xlabel('h: invagination distance/vesicle diameter')
           10     plt.savefig(fig_dir+'%s-vs-H_i.png' %(variable))
      
      Cell In [12], line 9, in plot_simulation(x, output, H_i, sigma, A_t, S_i, K, C_i, gamma_d, gamma_m, gamma_i, V_d, P, label, variable_df)
            7 y = outputs[output]
            8 plt.plot(x,y,label=label)
      ----> 9 plt.ylabel(variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1073, in _LocationIndexer.__getitem__(self, key)
         1070 axis = self.axis or 0
         1072 maybe_callable = com.apply_if_callable(key, self.obj)
      -> 1073 return self._getitem_axis(maybe_callable, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1312, in _LocIndexer._getitem_axis(self, key, axis)
         1310 # fall thru to straight lookup
         1311 self._validate_key(key, axis)
      -> 1312 return self._get_label(key, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1260, in _LocIndexer._get_label(self, label, axis)
         1258 def _get_label(self, label, axis: int):
         1259     # GH#5567 this will fail if the label is not present in the axis.
      -> 1260     return self.obj.xs(label, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/generic.py:4056, in NDFrame.xs(self, key, axis, level, drop_level)
         4054             new_index = index[loc]
         4055 else:
      -> 4056     loc = index.get_loc(key)
         4058     if isinstance(loc, np.ndarray):
         4059         if loc.dtype == np.bool_:
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexes/base.py:3805, in Index.get_loc(self, key, method, tolerance)
         3803     return self._engine.get_loc(casted_key)
         3804 except KeyError as err:
      -> 3805     raise KeyError(key) from err
         3806 except TypeError:
         3807     # If we have a listlike key, _check_indexing_error will raise
         3808     #  InvalidIndexError. Otherwise we fall through and re-raise
         3809     #  the TypeError.
         3810     self._check_indexing_error(key)
      
      KeyError: 'bending_energy_flipped'
      

      7959acb7198747215b87923bc4e915a713e285f5.png 600e27fff0050ff25e384a72c2918a3df9020746.png e4ee70b62b6cf5e183edc311b5d66452d9b178b3.png cf90324bc11aca34b40ec9377990dc45a2096708.png 3d8f5c9ea6d341d6c9370be3df24b175b924ed79.png c884cd3f4ce6c600fe2aaaeccc2b70af5c56173a.png 18711b65f3406c00f63c5929dc1a1c37f5fd5cc1.png d97eb35c03d5242069514d5d0ce9c507410a3143.png eb3dffa9a56b13b0f28f63e613bc26527a0c44b0.png 2a9847e2cce5327ddfeaea7e50875861364ab840.png 0bf9ae2988f6e83f580178f4cb126d0516bd88f9.png 8d8ddb1068bd7e69c6d139f79321222bd54ecaba.png 1e30ae3907710e74b76128b8f9a3cd34613aa9e3.png 5ad71de480c3dc0d19cb9d827176351e0b7babf3.png 28bd7517f43d52d58c90c7e104544a08fc4bebec.png

    2. as subplots of one figure
      outputs = llps_model()
      outputs.pop('theta_d')
      variables = list(outputs.keys())
      nplots = len(variables)
      width = 3
      height = nplots//3+1
      plt.figure(figsize=(width*5,height*4))
      
      for i, variable in enumerate(variables):
          plt.subplot(height, width, i+1)
          plot_simulation(output=variable)
          plt.xlabel('h: invagination distance/vesicle diameter')
      
      plt.tight_layout()
      plt.savefig(fig_dir+'all_outputs.png')
      
      ---------------------------------------------------------------------------
      KeyError                                  Traceback (most recent call last)
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexes/base.py:3803, in Index.get_loc(self, key, method, tolerance)
         3802 try:
      -> 3803     return self._engine.get_loc(casted_key)
         3804 except KeyError as err:
      
      File ~/.local/lib/python3.10/site-packages/pandas/_libs/index.pyx:138, in pandas._libs.index.IndexEngine.get_loc()
      
      File ~/.local/lib/python3.10/site-packages/pandas/_libs/index.pyx:165, in pandas._libs.index.IndexEngine.get_loc()
      
      File pandas/_libs/hashtable_class_helper.pxi:5745, in pandas._libs.hashtable.PyObjectHashTable.get_item()
      
      File pandas/_libs/hashtable_class_helper.pxi:5753, in pandas._libs.hashtable.PyObjectHashTable.get_item()
      
      KeyError: 'bending_energy_flipped'
      
      The above exception was the direct cause of the following exception:
      
      KeyError                                  Traceback (most recent call last)
      Cell In [43], line 11
            9 for i, variable in enumerate(variables):
           10     plt.subplot(height, width, i+1)
      ---> 11     plot_simulation(output=variable)
           12     plt.xlabel('h: invagination distance/vesicle diameter')
           14 plt.tight_layout()
      
      Cell In [12], line 9, in plot_simulation(x, output, H_i, sigma, A_t, S_i, K, C_i, gamma_d, gamma_m, gamma_i, V_d, P, label, variable_df)
            7 y = outputs[output]
            8 plt.plot(x,y,label=label)
      ----> 9 plt.ylabel(variable_df.loc[output]['description']+' ('+variable_df.loc[output]['units']+')')
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1073, in _LocationIndexer.__getitem__(self, key)
         1070 axis = self.axis or 0
         1072 maybe_callable = com.apply_if_callable(key, self.obj)
      -> 1073 return self._getitem_axis(maybe_callable, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1312, in _LocIndexer._getitem_axis(self, key, axis)
         1310 # fall thru to straight lookup
         1311 self._validate_key(key, axis)
      -> 1312 return self._get_label(key, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexing.py:1260, in _LocIndexer._get_label(self, label, axis)
         1258 def _get_label(self, label, axis: int):
         1259     # GH#5567 this will fail if the label is not present in the axis.
      -> 1260     return self.obj.xs(label, axis=axis)
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/generic.py:4056, in NDFrame.xs(self, key, axis, level, drop_level)
         4054             new_index = index[loc]
         4055 else:
      -> 4056     loc = index.get_loc(key)
         4058     if isinstance(loc, np.ndarray):
         4059         if loc.dtype == np.bool_:
      
      File ~/.local/lib/python3.10/site-packages/pandas/core/indexes/base.py:3805, in Index.get_loc(self, key, method, tolerance)
         3803     return self._engine.get_loc(casted_key)
         3804 except KeyError as err:
      -> 3805     raise KeyError(key) from err
         3806 except TypeError:
         3807     # If we have a listlike key, _check_indexing_error will raise
         3808     #  InvalidIndexError. Otherwise we fall through and re-raise
         3809     #  the TypeError.
         3810     self._check_indexing_error(key)
      
      KeyError: 'bending_energy_flipped'
      
      2203927067b93b2f1722c1864e9b15257a9ca1bb.png
  2. change in all outputs vs. invagination height for default parameters
    1. as individual plots
      outputs = llps_model()
      outputs.pop('theta_d')
      variables = outputs.keys()
      
      for variable in variables:
          plt.figure()
          plot_simulation_delta(output=variable)
          plt.tight_layout()
          plt.savefig(fig_dir+'%s-vs-H_i_delta.png' %(variable))
      

      7c4e2e47daff7d306cdf551f34dbcf9868a5ff3d.png 4b6697a24cb2dcd944390e68336e73bc662cf007.png 11b39b9fb19678cbeb7edf4416b604f62f6e92b4.png faa747541f5fe58cdae9eb8719b859838d358248.png bd73778a5ecded39f177dbb8e875e606e06a2dec.png 4b6697a24cb2dcd944390e68336e73bc662cf007.png d2e40413dd49c3476c891eb1b7a95a96edbee21b.png c081ae98e8e37c52701c3cd733a45928520b959b.png 5f82338248fc5b91a806502fffbd3f99af4f3d58.png cd58e0970412286596d706dfbbd7819216d7d2ee.png 2dda711a51a3e11b32a3b6d5f579bb5ed4bb005d.png 6d7b6d820619a81b8f26f5bb65e4e6561f4945fd.png 3e939816e92f8c16f4aceebe0be95cbc6469ca4a.png ef0af2a1669ad2f8260da613d5287aeb8a4c8072.png 5133f7b2ac66b87f7cdc48650905334741e87b70.png 0277a8d2f5285adb7488c714cad3590a78db5243.png 751e32f47c3756522b9af2225746363bec08600d.png 5375cb9b9242a08080972df42608b46ff7a26eee.png fdc8c78ac0f8c0b28099b190ac8cd2a664324874.png 3c5f25f93b4326c234cdafde725f704a9ae02bd9.png

    2. as subplots of one figure
      outputs = llps_model()
      outputs.pop('theta_d')
      variables = outputs.keys()
      nplots = len(variables)
      width = 3
      height = nplots//3+1
      plt.figure(figsize=(width*5,height*4))
      
      for i, variable in enumerate(variables):
          plt.subplot(height, width, i+1)
          plot_simulation_delta(output=variable)
      
      plt.tight_layout()
      plt.savefig(fig_dir+'all_outputs_delta.png')
      
      d2db7f07cc419da83f933f7e323ca28c863c779b.png
  3. change in surface areas
    variables = ['S_o', 'S_d', 'S_m']
    
    for variable in variables:
        plot_simulation_delta(output=variable, label=variable)
        plt.ylabel('∆ surface area')
        plt.legend()
        plt.tight_layout()
        plt.savefig(fig_dir+'areas-vs-H_i_delta.png')
    
    d0d49eb77b0368525e2d6cb1b6e45a61917cbc15.png
  4. each individual energy term output
    1. set temporary variables
      P_temp = 1e-3
      gamma_d_temp = gamma_d*10
      gamma_m_temp = gamma_m*10
      K_temp = 160
      

      :results: nil:END:

    2. normal droplet
      show_geometry(H_i=100)
      plt.savefig(fig_dir+'normal_droplet.png')
      
      2d3471c918354444e5dbfed3090b9a23fc8c084c.png
      1. with appropriate signs

        droplet only:

        variables = ['surface_tension', 'wetting_energy', 'droplet_energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'droplet_energy_contributions.png')
        plt.savefig(fig_dir+'droplet_energy_contributions.svg')
        
        a722597b21a9c4ef380799b145730cce9ccefe4a.png

        all contributions:

        variables = ['surface_tension', 'wetting_energy', 'droplet_energy',
                     'membrane_tension', 'bending_energy', 'turgor_pressure', 'F']
        colors = ['blue', 'green', 'black',
                  'red', 'purple', 'brown', 'orange']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=50, color=color)
        
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'energy_contributions.png')
        plt.savefig(fig_dir+'energy_contributions.svg')
        
        ccf25d920a31b9ca4e82423cf44185b73b83c59e.png
      2. flipped sign of negative energy contributions
        variables = ['membrane_tension', 'bending_energy_flipped','surface_tension',
                      'wetting_energy_flipped', 'turgor_pressure']
        
        plt.figure(figsize=(10,10))
        
        for variable in variables:
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'energy_contributions_flipped.png')
        
        d4fadf63b7f28f6cadec2250393c166db83137ef.png
    3. huge droplet
      show_geometry(H_i=100, V_d=V_d*10**4,
                    left=-3600, right=3600, bottom=-2300)
      plt.savefig(fig_dir+'bigdrop.png')
      
      6537c84d42334805fec7c558a6a2a6b65346a259.png
      1. with appropriate signs
        variables = ['surface_tension', 'wetting_energy', 'droplet_energy']
        colors = ['blue', 'green', 'black']
        
        plt.figure(figsize=(6,6))
        
        for variable, color in zip(variables, colors):
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, V_d=V_d*10**4, color=color)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'droplet_energy_contributions_bigdrop.png')
        plt.savefig(fig_dir+'droplet_energy_contributions_bigdrop.svg')
        
        43f8ed9f84d96dcd8f6a7b470c6ee9ace1a25127.png
      2. flipped sign of negative energy contributions
        variables = ['membrane_tension', 'bending_energy_flipped','surface_tension',
                      'wetting_energy_flipped', 'turgor_pressure']
        
        plt.figure(figsize=(10,10))
        
        for variable in variables:
            plot_simulation_delta(output=variable, label=variable, P=P_temp, gamma_d=gamma_d_temp,
                                  gamma_m=gamma_m_temp, K=K_temp, V_d=V_d*10**4)
        plt.legend()
        #plt.ylim(-100, 100)
        plt.ylabel('∆ energy (pN * nm)')
        #plt.title(title)
        plt.tight_layout()
        plt.savefig(fig_dir+'energy_contributions_bigdrop_flipped.png')
        
        6d30e3a0070ec0b109183d4eeb7b92abcd332144.png
  5. Plot F for a series of values of gammam

    What are the parameter values? I assume that gammam is negative. Plot F for a series of values of gammam. Based on the plots for Sm, Sd, and S0, F should end up with a minimum at H/R=1 if gammam is large enough. It might also happen if gammad is small enough. It's hard to say because I can't easily compare the magnitude of changes in Sd, Si, and S0.

    1. absolute
      gamma_m_range = [0,50,100,150,199]
      
      for i in gamma_m_range:
          plot_simulation(gamma_m=i, label=i/gamma_d)
      
      plt.legend(title='gamma_m/gamma_d')
      plt.xlabel('invagination distance/vesicle diameter')
      plt.tight_layout()
      plt.savefig(fig_dir+'f_scan_gamma_m.png')
      
      /tmp/ipykernel_195576/574978439.py:6: RuntimeWarning: invalid value encountered in arccos
        theta_d = np.arccos(gamma_m/gamma_d) #equation 8
      
      b5c67e459db1961a1ac42be598b9f812e20b5cd2.png
    2. normalized
      gamma_m_range = [0,50,100,150,199]
      
      for i in gamma_m_range:
          plot_simulation_normalized(gamma_m=i, label=i/gamma_d)
      
      plt.legend(title='gamma_m/gamma_d')
      plt.xlabel('invagination distance/vesicle diameter')
      plt.tight_layout()
      plt.savefig(fig_dir+'f_scan_gamma_m.png')
      
      /tmp/ipykernel_195576/574978439.py:6: RuntimeWarning: invalid value encountered in arccos
        theta_d = np.arccos(gamma_m/gamma_d) #equation 8
      
      ba9679aa27cc8ccb109502ceee26f62ab0e3ffb2.png

      :end:

    3. delta
      gamma_m_range = np.array( [0.01,0.25,0.5,0.75] )*gamma_d
      
      for i in gamma_m_range:
          plot_simulation_delta(gamma_m=i, label=i/gamma_d)
      
      plt.legend(title='$\gamma_{m}/\gamma_{d}$')
      plt.tight_layout()
      plt.savefig(fig_dir+'f_scan_gamma_m_delta.png')
      
      57d4efde7c67d058f2e8b3a6df01e3a25041b29a.png
    4. delta and geometry side by side
      gamma_m_range = [0,50,100,150,180]
      
      for i in gamma_m_range:
          plot_simulation_geometry(gamma_m=i, label=i/gamma_d)
          plt.legend(title='$\gamma_{m}/\gamma_{d}$')
          plt.tight_layout()
          plt.savefig(fig_dir+'f_scan_gamma_m_%s.png' %(i))
      
      /tmp/ipykernel_195576/574978439.py:6: RuntimeWarning: invalid value encountered in arccos
        theta_d = np.arccos(gamma_m/gamma_d) #equation 8
      /home/maxferrin/miniconda3/envs/llps_model_analysis/lib/python3.10/site-packages/matplotlib/axes/_base.py:2480: UserWarning: Warning: converting a masked element to nan.
        xys = np.asarray(xys)
      

      0b4e1e340996085e06113afc9baae0ec7c367b42.png 5c87a9bf8033b4ecb0b2e64b312beb29e43f6515.png 7b656cc1cdc1ec41f7f9f4f018ce9e27c8771db6.png 426b73223917a8f8e24977be236d40e9323b9300.png 36df76ae297ce9802abfd20abaac9264a7dbfd81.png

  6. F over series of turgor pressure
    pressure_range = np.linspace(10**-2,1,5)
    
    for i in pressure_range:
        plot_simulation_delta(P=i, label=i)
    
    plt.legend(title='turgor pressure ($N/nm^{2}$)')
    plt.tight_layout()
    plt.savefig(fig_dir+'f_scan_turgor_delta.png')
    
    7c48277dacb237516228350504e0932d7cafcc87.png
  7. spontaneity analysis

    there's no solution for kinetic \(\sigma\)

    1. DONE \(V_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(7,9,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(V_d=x_range)
          kinetic_boundary = kinetic_boundaryy(V_d=x_range)
          both_boundary_test(x_param='V_d', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.9,1,1.1])
          plt.savefig(fig_dir+y_param+'V_d_boundary_check.png')
      

      543ac80568d7e04edfb71d20fd0720da69372f9d.png dfd95ca95ba7a46a85c98464a3b30d0af800497b.png 7e6f6ee0c723be8c9cf077817d05bf557f6d3d37.png

      show_geometry(H_i=100, V_d=10**(7))
      plt.savefig(fig_dir+'10e7_geometry.png')
      
      0194c353098ccd34b89dfc2b7f5a07b19ce24989.png
      show_geometry(H_i=100, V_d=10**(9), left=-1200, right=1200, bottom=-600)
      plt.savefig(fig_dir+'10e9_geometry.png')
      
      81999fed97c82d0878e965dc4a4717bcb83864db.png
    2. DONE \(\gamma_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-4,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_d=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_d=x_range)
          both_boundary_test(x_param='gamma_d', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.9,1,1.1])
          plt.savefig(fig_dir+y_param+'gamma_d_boundary_check.png')
      

      1c5171953243969544a193cbabd0fbf074a10282.png c0e2ddd0ce7df6a32043f779bfbd9d09c87c844e.png 3fee03a2a46f560870850b2f20a0fd89716606a9.png

    3. DONE \(\gamma_{m}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-4,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_m=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_m=x_range)
          both_boundary_test(x_param='gamma_m', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*2, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'gamma_m_boundary_check.png')
      

      06e3aae43985001e6cdcbca291a145f76f7c7064.png 4fdf88ffa9b8b5142867b9836692b8decb089a60.png bbedad9e0ef40f736f0309e7e6e66b0767101583.png

    4. TODO \(\gamma_{m}\) /\(\gamma_{d}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    #(thermo_spont_sigma,kinetic_spont_sigma, 'sigma'),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-2,0.0001,base=10)*gamma_d
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(gamma_m=x_range)
          kinetic_boundary = kinetic_boundaryy(gamma_m=x_range)
          both_boundary_test(x_param='gamma_m', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/2, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'gamma_m_relative_boundary_check.png')
      

      22e6df76388400797c3c170e9b5881786d8de935.png b9612dd8906b097f2834abb0d65b5821118c69ce.png 3c00992cdbfd2a83f4c19f7ef297a68835970678.png

    5. DONE \(\sigma\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', -2),
                    (thermo_spont_P,kinetic_spont_P, 'P', 2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', -2)]
      x_range = np.logspace(-5,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(sigma=x_range)
          kinetic_boundary = kinetic_boundaryy(sigma=x_range)
          both_boundary_test(x_param='sigma', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/5, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'sigma_boundary_check.png')
      

      cbb28c5779836d223f04ec7bcafe0818bc4d637a.png 46808b0992afe4806d800ac432cd0ab72d982d63.png 90f72c29e56e420ad8edd069dede3c2fea81c384.png

    6. DONE \(S_{i}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(3.5,4.5,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(S_i=x_range)
          kinetic_boundary = kinetic_boundaryy(S_i=x_range)
          both_boundary_test(x_param='S_i', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*50, factors=[0.999,1,1.001])
          plt.savefig(fig_dir+y_param+'S_i_boundary_check.png')
      

      766d95fc08012c0180527a485b42a44a3768eb0a.png f1c90e32615ae1b3bd89ee522df0ce3092706962.png 5b48516fbf16ca43d24f4eb13305d3747c718c8b.png

    7. DONE \(P\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    #(thermo_spont_P,kinetic_spont_P, 'P', 2),
                    (thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2)]
      x_range = np.logspace(-3,0,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(P=x_range)
          kinetic_boundary = kinetic_boundaryy(P=x_range)
          both_boundary_test(x_param='P', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset*150, factors=[0.999,1,1.001])
          plt.savefig(fig_dir+y_param+'P_boundary_check.png')
      

      fe6f00b4583ce1f366f083c1b8bb998c671212f0.png 955ba4e0d0369cf9f356a3f397c8b22aae0ffb7a.png

    8. DONE \(C_{i}\) x-axis
      boundaries = [(thermo_spont_K,kinetic_spont_K, 'K', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2)]
      x_range = np.logspace(-4,-1,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(C_i=x_range)
          kinetic_boundary = kinetic_boundaryy(C_i=x_range)
          both_boundary_test(x_param='C_i', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/5, factors=[0.8,1,1.2])
          plt.savefig(fig_dir+y_param+'C_i_boundary_check.png')
      

      cf3dff2f28859a18dc0840d2defd0483142f702e.png c861a6b00b033dd106709caf61eb444144ad9725.png

    9. DONE \(K\) x-axis
      boundaries = [(thermo_spont_C_i,kinetic_spont_C_i, 'C_i', 2),
                    (thermo_spont_P,kinetic_spont_P, 'P', -2)]
      x_range = np.logspace(1,4,base=10)
      for thermo_boundaryy, kinetic_boundaryy, y_param, y_offset in boundaries:
          thermo_boundary = thermo_boundaryy(K=x_range)
          kinetic_boundary = kinetic_boundaryy(K=x_range)
          both_boundary_test(x_param='K', y_param=y_param, x_range=x_range,
                          thermo_boundary=thermo_boundary, kinetic_boundary=kinetic_boundary,
                          y_offset=y_offset/1, factors=[0.9,1,1.1])
          plt.savefig(fig_dir+y_param+'K_boundary_check.png')
      

      2ee978fd97452564c699b6151762d86491451714.png fb55d20816eea2679944def168db4db6c0bba592.png

  8. plots for 20200915 lab meeting drafting
    1. TODO each individual energy term output
      variables = ['membrane_tension',
                    'bending_energy','surface_tension',
                    'wetting_energy','turgor_pressure']
      
      for variable in variables:
          plt.figure()
          plot_simulation_delta(output=variable)
          plt.tight_layout()
          plt.savefig(fig_dir+'%s-vs-H_i_delta.png' %(variable))
      

      3e939816e92f8c16f4aceebe0be95cbc6469ca4a.png ef0af2a1669ad2f8260da613d5287aeb8a4c8072.png 0277a8d2f5285adb7488c714cad3590a78db5243.png 751e32f47c3756522b9af2225746363bec08600d.png 3c5f25f93b4326c234cdafde725f704a9ae02bd9.png

    2. membrane components without droplet
      1. effect of membrane tension alone
        plot_simulation_geometry(V_d=0, K=0, gamma_m=0, gamma_i=0, P=0)
        
        9ea78939015affe7240ce9cb3b53c6217484983d.png
      2. effect of bending energy alone
        plot_simulation_geometry(K=K, V_d=0, sigma=0, gamma_m=0, gamma_i=0, P=0)
        
        cb5ea66c30b957fa8608bddd344b6affd3141b1c.png
        1. TODO fix this

          why is this different from plotting membrane bending energy term alone?

          R_i_range = np.linspace(0.1,np.pi,100)
          bending_energy = (K/2)*S_i*((2/R_i_range)-C_i)**2
          plt.plot(R_i_range, bending_energy)
          
          <matplotlib.lines.Line2D at 0x7ff9714e7ac0>
          e68abfd1a0ec650df791468df1e56efc7eecf549.png
    3. droplet effects
      1. surface tension
        plot_simulation_geometry(sigma=0, K=0, gamma_m=0, gamma_i=0, P=0)
        
        53e82532d5d2f8347b71a50ebaefa884bb7c4a3c.png
      2. surface tension + wetting energy

        can't have wetting energy alone

        plot_simulation_geometry(sigma=0, K=0, P=0)
        
        ac28dcb5037d3e31bd1c10f1b142f79ec574c0d8.png
    4. TODO droplet+membrane effects
    5. TODO turgor pressure effects

Author: Max Ferrin

Created: 2023-01-23 Mon 15:25

Validate